© ESO, 2017

1. Introduction

Observing the circum-galactic medium (CGM) represents an important challenge for understanding how galaxies form and evolve. Galaxy evolution is driven primarily by the flows of gas that surround galaxies. Moreover, the CGM contains a large amount of the baryonic matter in galaxies and as such observations of this gas provide crucial information. A powerful tracer of this gas is Lyman alpha (Lyα) emission, which allows circum-galactic gas to be observed around high-redshift galaxies as a Lyα halo. A number of physical mechanisms can contribute to spatially extended Lyα emission, including fluorescence, cooling radiation, or the scattering of Lyα photons produced in star-forming HII regions (Gould & Weinberg 1996; Katz et al. 1996; Haiman et al. 2000; Haiman & Rees 2001; Cantalupo et al. 2005; Dijkstra et al. 2006; Kollmeier et al. 2010; Barnes & Haehnelt 2010; Lake et al. 2015).

Such extended Lyα emission has been detected around nearby galaxies (e.g. Kunth et al. 2003; Hayes et al. 2005; Hayes 2015). By selecting 14 nearby galaxies that cover the same range of far-UV luminosities as high-z galaxies, the Lyα Reference Sample (LARS; Östlin et al. 2014) collaboration constructed a sample that is comparable to high-z samples. Most of their galaxies show Lyα emission that is more extended than both the stellar UV continuum and the Hα emission showing the rich gas content of the CGM (Hayes et al. 2013, 2014; Herenz et al. 2016).

At high redshift, the mapping of the extended Lyα haloes around galaxies (non-AGN) is however a lot more difficult because of sensitivity and resolution limitations. Detections of extended Lyman alpha emission at high redshift have been obtained in the past. While some large Lyα blobs have been observed (e.g. Steidel et al. 2000; Matsuda et al. 2004, 2011), most of these studies were forced to employ stacking analyses because of sensitivity limitations. The first tentative detections of Lyα haloes around normal star-forming galaxies emitting Lyα emission using narrowband (NB) imaging methods were reported by Møller & Warren (1998) and Fynbo et al. (2001). Later, Hayashino et al. (2004) observed 22 Lyman break galaxies (LBG) and detected extended Lyα emission by stacking the NB images. These authors were followed six years later by Ono et al. (2010) who detected Lyα haloes in their composite NB images of 401 Lyα emitters (LAEs) at z = 5.7 and 207 at z = 6.6. Matsuda et al. (2012) and Momose et al. (2014) significantly increased the size of LAEs samples used by stacking ≈2000 and ≈4500 LAEs at redshift z ≃ 3 and 2.2 ≤ z ≤ 6.6, respectively. Momose et al. (2014) found typical Lyα halo exponential scale lengths of 5–10 physical kpc. Matsuda et al. (2012) found that Lyα halo sizes are dependent on environmental density; these halo sizes extend from 9 to up to 30 physical kpc towards overdense regions. More recently, Xue et al. (2017) studied ≈1500 galaxies in two overdense regions at z ≈ 3 and 4. Using stacking methods these authors reported Lyα halo exponential scale lengths of 5–6 physical kpc and found that Lyα halo sizes correlate with the UV continuum and Lyα luminosities, but not with overdensity. Steidel et al. (2011) stacked 92 brighter (R_AB ≃ 24.5) and more massive LBGs at z = 2.3−3, finding large Lyα extents of ≈80 physical kpc beyond the mean UV continuum size at a surface brightness level of ~10^-19 erg s^-1 cm^-2 arcsec^-2. Put together, all these studies showed that Lyman alpha emission is on average more spatially extended than the UV stellar continuum emission from galaxy counterparts.

Meanwhile, other studies have found conflicting results. Feldmeier et al. (2013) argued that the observed extended emission is artificially created by an underestimation of the stacking procedure systematics. After carrying out an error budget analysis, they did not find evidence for significant extended Lyα emission. Bond et al. (2010) also reported compact Lyα emission in their stack of eight star-forming galaxies at z = 3.

Over a similar period and using a different approach, Rauch et al. (2008) performed an ultra-deep (92h) long-slit observation and identified 27 faint LAEs (few ×10^-18$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$) at redshift 2.67 <z< 3.75. This observation enabled the individual detections of extended Lyα emission along the slit for most of their objects although with large uncertainties owing to slit losses and the high errors on the continuum size measurements. Some other detections of extended Lyα emission around high-redshift star-forming galaxies were obtained using the magnification power of gravitational lensing (e.g. Swinbank et al. 2007; Patrício et al. 2016).

Recently, a significant step forward has been taken thanks to the substantial increase in sensitivity provided by the Multi-Unit Spectroscopic Explorer (MUSE) at the ESO-VLT (Bacon et al. 2010). Wisotzki et al. (2016; hereafter W16) reported the detection of 21 Lyα haloes around relatively continuum-faint (m_AB ≳ 27) star-forming galaxies at redshift 3 <z< 6 within the Hubble Deep Field South (HDFS) observed with MUSE. Their data reach an unprecedented limiting surface brightness (SB) of ~10^-19ergs^-1cm^-2arcsec^-2 (1σ) enabling the study of the CGM on a galaxy-by-galaxy basis. The Lyα haloes from the W16 study have exponential scale lengths ranging from ≈1 kpc to ≈7 kpc and appear to be on average 10 times larger than their corresponding UV galaxy sizes. These new observational data also enable the direct comparison of the Lyα halo properties with the stellar properties of the host galaxies and the investigation of the origin of the Lyα haloes. This pioneering study was however limited to a small sample and therefore the results need to be confirmed with better statistics.

Here, we extend the W16 LAE sample by one order of magnitude using the Hubble Ultra Deep Field (UDF) data obtained with MUSE (Bacon et al. 2015). The significant effort on the data reduction of this data set improves the limiting SB sensitivity by one order of magnitude over previous narrowband studies. First, we follow a similar approach as W16 to quantitatively characterize the spatial extent of the Lyα emission around high-redshift galaxies in the UDF (− 15 ≥ M_UV ≥ −22). We then analyse the sizes and Lyα luminosities of our Lyα haloes as a function of the UV properties of their HST counterparts and compare our results to W16. In addition to its spatial distribution, the Lyα line profile encodes crucial information that can help shed light on the origin of the Lyα emission and constrain the gas opacity and kinematics (Haiman et al. 2000; Dijkstra et al. 2006; Verhamme et al. 2006; Kollmeier et al. 2010; Gronke & Dijkstra 2016). Taking advantage of the spectral information of MUSE data cubes, we also investigate how Lyα emission relates to various line properties, such as the line width and equivalent width.

The paper is organized as follows: we describe our data and our sample construction in Sect. 2. Section 3 presents our procedure for the extraction of the images and construction of radial SB profiles needed for the detection of extended Lyα emission. Section 4 explains the Lyα spatial distribution modelled that we use to determine the characteristics of the Lyα haloes that are presented in Sect. 5. Section 5 also includes the analysis of the Lyα line profile. In Sect. 6 we investigate the relation between the Lyα haloes and their host galaxies. Finally, we discuss our results in Sect. 7 and present our summary and conclusions in Sect. 8. Appendix A gives a comparison of the Lyα haloes detected around galaxies, which are both in the deep udf-10 data cube and in the shallower mosaic data cube.

For this paper, we use AB magnitudes, physical distances, and assume a ΛCDM cosmology with Ω_m = 0.3, Ω_Λ = 0.7, and H₀ = 70 km s^-1 Mpc^-1.

2. Data and sample definition

2.1. Observations and data reduction

The UDF data were taken using the MUSE instrument between September 2014 and February 2016 under the MUSE consortium Guarantee Time Observations. A number of 1′ × 1′ pointings (corresponding to the MUSE field of view) were completed at two levels of depth. The medium-deep data consist of a mosaic of 9 deep, 10 h pointings denoted udf-0[1-9]. The ultra-deep data, denoted udf-10,  consist of a single ≈20 h pointing that overlaps with the mosaic reaching a total of 30 h depth. During the observations the sky was clear with good seeing conditions (full width at half maximum (FWHM) of at 7750 Å). More details about the data acquisition can be found in Bacon et al. (2017; hereafter B17).

The data reduction of both the udf-10 and mosaic data cubes is described in B17. The two resulting data cubes contain 323 × 322 and 945 × 947 spatial pixels for the udf-10 and mosaic field, respectively. The number of spectra match the number of spatial pixels in each data cube with a wavelength range of 4750 Å to 9350 Å (3681 spectral pixels) with medium spectral resolution R ~ 3000. The spatial sampling is per spaxel and the spectral sampling is 1.25 Å per pixel. The data cubes also contain the estimated variance for each pixel. The data reach a limiting SB sensitivity (1σ) of 2.8 and 5.5 × 10^-20 erg s^-1 cm^-2Å^-1 arcsec^-2 for an aperture of 1′′ × 1′′ in the 7000–8500 Å range for the udf-10 and mosaic data cubes, respectively (see B17 for more details).

Based on these reduced data cubes we constructed two catalogues corresponding to each data cube. The source detection and extraction were performed using HST priors, imposing a magnitude cut at 27 in the F775W band for the mosaic field only, and the ORIGIN (Mary et al., in prep.) detection software. A complete description of the strategy used for the catalogue construction can be found in Inami et al. (2017; hereafter I17). The ORIGIN software (see B17 for technical details) is designed to detect emission lines in 3D data sets. This software enables the discovery of a large number of LAEs that are barely seen or even undetectable in the HST images. Photometric magnitudes for these new objects were calculated following the method described in B17.

2.2. Lyman alpha emitters sample

Our parent sample was constructed from UDF catalogues (see I17) and according to the following criteria:

• 1.

  We selected the LAEs (“TYPE = 6” in the catalogues) with areliable redshift (“CONFID = 2 and 3”). This yields a sample of155 and 620 objects for the udf-10 and mosaic,respectively.

• 2.

  Our primary objective being the study of individual galaxies, we removed galaxies in pairs closer than 50 kpc in projected transverse separation and with velocity differences of less than 1000 km s^-1 , which was estimated using the peak of the Lyα line or the red peak if the line was double peaked. We found 28 and 64 such objects in the udf-10 and mosaic, respectively. The study of the Lyα haloes of such LAE pairs will be part of another study. The analysis of merger rates from the MUSE UDF data is detailed in Ventou et al. (2017).

• 3.

  We also excluded 20 and 25 objects that are closer than 3′′ and 4′′ to the edges of the udf-10 and mosaic data cubes, respectively. This is necessary to ensure we can analyse extended Lyα emission over a large spatial window for our entire sample. Objects from the udf-10 data cube are allowed to be closer to the edges because of the higher quality of the data given that the udf-10 data cube is combined with the wider mosaic data cube.

• 4.

  Among the remaining objects, we manually removed 7 and 29 objects in the udf-10 and mosaic fields, respectively, which are contaminated by emission lines from foreground sources, skyline residuals, or by continuum remnants visible in the NB image (see Sect. 3.1.1 for the continuum subtraction method).

• 5.

  Finally, following the procedure described in Sect. 3.1.1, we created NB images around the Lyα emission line and imposed a minimal signal-to-noise ratio (S/N) of 6 in a fixed and large aperture set using a curve of growth (CoG) method (see Sect. 5.3.2). The S/N is defined as the Lyα flux divided by the standard deviation from the data cube. This cut is motivated by our detection limit estimation described in Sect. 4.3.1. It eliminates 43 LAEs in the udf-10 field and 282 in the mosaic field. The S/N cut introduces a selection bias towards brighter haloes. This bias is noticeable in the lower panel of Fig. 1, where the total Lyα flux distribution before and after the S/N cut is shown.

In total 26 galaxies are in both udf-10 and mosaic fields. For these objects, we only show results from the udf-10 data cube because of the higher S/N; a comparison of the results from the two data cubes is given in Appendix A. Our final sample consists of 252 galaxies: 57 in the udf-10 field and 195 that are uniquely in the mosaic field. The sample spans a redshift range from 2.93 to 6.04 and a total Lyα flux ranging from ≈ 1.6 × 10^-18$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$ to ≈ 1.1 × 10^-16$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$. Lyα fluxes are measured using a CoG method (see Sect. 5.3.2). The redshift and flux distributions of our sample (purple) and the sample without the S/N cut (grey) are shown in Fig. 1. The flux distribution shows the selection bias towards the brighter LAEs.

3. Detection of diffuse Lyα emission

Our detection of extended Lyα emission employs a circularly symmetric analysis, following a similar approach to W16. The method uses the radial SB profiles of both the Lyα and UV continuum emission. In this section we first describe the methods used to create the Lyα NB and UV continuum images, then we explain how the radial surface brightness profiles are constructed, and finally we present some LAE radial SB profiles as examples.

3.1. Image construction

3.1.1. Lyα narrowband images

We constructed a 10″ × 10″ Lyα NB image for each object from the MUSE data cube. We used a wide, fixed spatial aperture to ensure that we included all of the detectable Lyα emission around our galaxies. In order to remove the continuum, we first performed a spectral median filtering on the MUSE data cube in a wide spectral window of 200 spectral pixels, in effect removing any emission lines (see Herenz & Wisotzki 2017 for the validation of this continuum subtraction method on MUSE data cubes). A continuum-free data cube was then computed by subtracting the filtered data cube from the original. In some cases the continuum of very bright objects was not well subtracted. As specified in Sect. 2.2, the 21 affected objects were removed from the sample.

The spectral bandwidths of the Lyα NB images were defined to maximize the S/N in a fixed spatial aperture (radius of 2′′). Following this procedure, we obtained NB images with spectral bandwidths ranging from 2.5 to 20 Å (i.e. 2 and 16 MUSE pixels, respectively). The largest spectral bandwidths correspond to double-peaked lines (some examples can be seen in Figs. 2 and 3). The average spectral width is 6.25 Å.

3.1.2. Ultraviolet continuum images

We constructed UV continuum images for our sample using one of three different HST images of the UDF (Illingworth et al. 2013), depending on the redshift of the object. The F814W ACS/WFC, F105W WFC3/IR, and F125W WFC3/IR HST images are used for objects at z< 4, 4 ≤ z< 5, and z ≥ 5, respectively. We chose these filters because they are not contaminated by the Lyα emission or by intergalactic medium (IGM) absorption. These filters also probe UV continuum over similar rest-frame wavelength ranges, which are approximately 1400–2300 Å, 1500–2400 Å, and 1570–2300 Å for the F814W ACS/WFC, F105W WFC3/IR, and F125W WFC3/IR HST filters, respectively.

For each object in our sample, we constructed UV continuum images with the HST counterparts from the I17 catalogue. After masking the pixels outside the segmentation map for each HST counterpart, we resampled the masked HST images to MUSE resolution and convolved them with the MUSE PSF. The HST PSF is not taken into account here because, first, the HST PSF (FWHM of 0.09′′ for the F814W band and 0.19′′ for the F105W and F125W bands – Rafelski et al. 2015) is much smaller than the MUSE PSF (≈0.7′′) and, second, the constructed UV continuum images are only used to compare visually Lyα and UV spatial extents. Our UV continuum modelled based on the HST data (see Sect. 4.1) considers the HST PSF.

The method used to estimate the wavelength-dependent PSF of the udf-10 and mosaic data cubes is detailed in B17. It is best described as a two-dimensional Moffat distribution with a fixed beta parameter of 2.8 and a wavelength-dependent FWHM that we evaluate at the wavelength of each Lyα line.

Twenty-one of our objects are not in the Rafelski et al. (2015) HST catalogue and instead were discovered by ORIGIN; however, these 21 objects are visible in the HST image. Magnitudes and segmentation maps for these galaxies were calculated using NoiseChisel (Akhlaghi & Ichikawa 2015) and added to the MUSE UDF catalogues (see B17 and I17). Thirteen other Lyα emitters of our sample discovered by MUSE do not show any HST counterpart (e.g. object #6498 in Fig. 2). When comparing to the corresponding Lyα radial SB profiles, we treat these galaxies as point-like sources convolved with the MUSE PSF.

We could have constructed continuum images directly from the MUSE data cubes. However, this is only possible for the brightest objects as most of our objects have poor continuum S/N in the MUSE data cubes. In addition, source blending is important at the MUSE spatial resolution while at HST resolution most of our sources are well separated.

3.2. Surface brightness radial profiles

To visually compare the spatial extents of the UV and Lyα emission, we constructed radial SB profiles. We performed aperture photometry on the Lyα NB images and UV continuum images by averaging the flux in successive, concentric, one-pixel-wide annuli centred on the Lyα emission centroid. For the objects in our sample without HST counterparts, we compared Lyα radial SB profiles to the MUSE PSF radial SB profiles. The Lyα centroid was measured by fitting a simple 2D Gaussian to the Lyα NB image. In some cases, the centroid measured from the MUSE data is offset from the coordinates from the HST catalogue. The offsets are relatively small: less than 0.3″ for ≈95% of our sample (median value ≲0.1″). We therefore ignored these offsets when constructing SB profiles and assumed the UV and Lyα emission to be concentric. Errors on Lyα radial SB profiles were measured in each annulus using the estimated variance from the MUSE data cubes.

Figures 2 and 3 show a representative subsample of 14 objects from the mosaic and udf-10 fields. These objects were chosen to exhibit the diversity of the LAEs in terms of luminosity, line profile, and spatial extent. For each object we show the corresponding HST image we used for the study (see Sect. 3.1.2); the MUSE white light image, summed over the full MUSE data cube spectral range; the Lyα line, which is integrated in the HST counterparts mask convolved with the MUSE PSF (see the white contours on the white light image); the Lyα NB image (see Sect. 3.1.1); and the radial SB profiles. The UV continuum and PSF profiles have been re-scaled to the Lyα emission profile to aid the visual comparison.

Most of the objects show Lyα emission that is more spatially extended than the UV continuum. Some objects display a clear Lyα halo (e.g. objects #1185, #82, #1087, #53, and #6297) but for other objects the extended Lyα emission is not as obvious (objects #6498, #6534, or #218). Further analysis of the statistical significance of the detected Lyα haloes is presented in Sect. 5.1.

4. Lyα halo modelling

4.1. Two-dimensional two-component fits

Here, we describe how we characterize the spatial distribution of extended Lyα emission. Following W16, we fit Lyα emission with a two-dimensional, two-component exponential distribution using the Python/photutils package (Bradley et al. 2016). Specifically, we decomposed the observed 2D Lyα distribution into central and extended exponential components using the HST morphological information as prior. The W16 work demonstrated that this decomposition is appropriate for characterizing Lyα haloes for a similar LAE sample. Adopting the same approach allows us to compare to their results directly.

The modelling is performed in two distinct steps:

• 1.

  First, the UV continuum is fit with a circular, 2D exponentialdistribution. We directly fit the HST image, chosen dependingon the redshift of the object, (see Sect. 3.1.2)taking into account the corresponding PSF (Rafelskiet al. 2015, Table 1). The HSTcounterpart of a given object was isolated by masking its sur-roundings using the HST segmentation mask. This first fit yieldsthe continuum spatial scale length of each host galaxy. The Rafel-ski et al. (2015) HST segmentationmask was created by combining the detection maps of the object inseveral HST bands. It is therefore supposed to delimit the galaxyin a rather large area and thus include most of the UV flux. If the ob-ject is located in a crowded region, the mask can be smaller to allowthe separation of the sources. This is however very rare becauseof the high resolution of the HST images. Moreover, we find thatthe extent of the HST segmentation mask has a small impact on theresulting scale length. This is because the fit is mainly driven bythe central emission of the galaxy. Consequently, even if there arepotential faint UV counterparts surrounding the galaxy, the scalelengths are not drastically different.

• 2.

  Second, the Lyα NB image is fit by a sum of two circular, 2D exponential distributions, fixing the scale length of the first component to the continuum distribution value. Hence, the first component corresponds to central, core emission and the second to emission from an extended halo. The fit takes into account the MUSE PSF by convolving the model with the PSF and the variance of each pixel in the image.

We thus have three parameters in total to fit the Lyα distribution: the halo scale length and fluxes of both the Lyα core and halo.

Figure 4 shows the best-fit model radial SB profiles, decomposed into core emission (green line), extended halo emission (blue line), and with the total emission shown in red. We overplot the Lyα SB radial profiles (black dots with error bars). For most of our objects, the modelled radial SB profiles are a good representation of the observed profiles. The 2D, two-component decomposition model therefore appears to be a good description of the Lyα distribution around LAEs. The W16 authors also found this decomposition to be consistent with their observed data.

4.2. Error estimation

We estimated errors on the best-fit halo scale length measurements. First, we generated 100 realizations of each best-fit model Lyα image by combining the noise-free model image with realizations of the estimated noise. The noise was assumed to follow a normal distribution with the variance at each pixel set equal to the variance of the corresponding pixel in the MUSE data cube. Each realization of a given object was then fit and the final error on the halo scale length was given by the standard deviation across the recovered scale lengths.

To estimate the error on the core scale length (which is instead fit to the HST image) we followed a similar procedure using 100 empty regions of the HST image as artificial noise.

4.3. Detection limit

4.3.1. Signal-to-noise limit

To estimate the limitations of our 2D decomposition, we fit a range of simulated Lyα distributions combined with random realizations of the noise again using the variance from the MUSE data cubes. The variance used here was estimated around 6000 Å in a 6 Å spectral window corresponding to the median NB image spectral bandwidth of our sample. This allows us to assess the S/N needed to measure Lyα halo properties reliably from our observed sample.

We considered simulated Lyα distributions with a fixed core scale length, rs_cont, a fixed core flux, F_cont, a range of 5 halo scale lengths, rs_halo and a broad range of halo fluxes, F_halo. The fixed core values were set to the averages of our sample, F_cont = 4.0 × 10^-18$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$ and rs_cont = 0.06′′ (i.e. 0.3 MUSE pixel). We considered halo fluxes ranging from 1 × 10^-20 to 2 × 10^-18$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$ and a set of halo sizes, rs_halo = [ 0.2′′, 0.4′′, 0.6′′, 1.0′′, 1.4′′, 1.8′′, 2.2′′] (i.e. [1, 2, 3, 5, 7, 9, 11] MUSE pixels). For each model Lyα distribution we generated 100 noise realizations and assessed the success rate for reliably recovering the halo size.

We find that halo sizes are reliably recovered above a S/N ≈ 6. The S/N is measured inside an aperture corresponding to the CoG radius r_CoG, which represents the radius for which the averaged flux in a concentric 1-pixel annulus reaches the noise value (see Sect. 5.3.2). The smaller Lyα haloes are therefore less penalized by a S/N cut estimated within this aperture than if we had used a wide aperture that is identical for every object.

We also considered a range of core values (rs_cont, F_cont, not shown here) and find that this value of the S/N limit is still appropriate. This S/N cut was thus adopted for the sample construction (see Sect. 2.2).

4.3.2. Size and flux limit

In addition to being limited in sensitivity by S/N, we are also limited in our ability to measure the sizes of very compact Lyα haloes by the MUSE PSF. To estimate the Lyα halo scale length below which we cannot trust our measurements, we again ran our modelled routine on several model Lyα distributions (with artificial noise based on the variance data cube). For each model object, we incrementally decreased the Lyα halo scale length until we could no longer recover the input value. We find that the resulting scale length limit corresponds to one quarter of the MUSE PSF FWHM. In practice, the halo scale length limit is a function of wavelength due to the PSF dependence on wavelength and thus on redshift. As such, we calculated the limit separately for each object in our observed sample, yielding values ranging from 0.85 kpc to 1.48 kpc. If the best-fit Lyα scale length was below the scale length limit, we considered the limit value as an upper limit.

We also tested our ability to detect the faint Lyα haloes reliably. We performed the same exercise as our S/N limit procedure (see previous subsection) and find as expected that our halo flux limit increases with the halo scale length. This is because the surface brightness shrinks as the total Lyα flux is preserved. Fixing the core component to the averages of our sample, we deduced that the halo flux limit corresponds to 9 × 10^-19$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$ and 5 × 10^-19$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$ times the halo scale length for the mosaic and udf-10, respectively.

Lyα halo measurements and fitting results are given in Table B.1.

5. Lyα halo characteristics

In this section, we present the characteristics of our sample of Lyα haloes in terms of sizes, fluxes, spatial, and spectral shapes.

5.1. Statistical significance of the detected Lyα haloes

Of the 252 galaxies, 87 objects either have Lyα halo fluxes (62 objects), scale lengths (19 objects), or both (6 objects) below the detection limits (see Sect. 4.3.2). For the 19 objects with only the Lyα halo scale length below the detection limit, we use this value as an upper bound. The objects with Lyα halo fluxes below the detection limit are ignored in the rest of this section. This leaves us with a sample of 184 LAEs at this stage.

Following W16, in order to evaluate the statistical significance of the detected extended Lyα emission, we calculated the probability p₀ of the two scale lengths (galaxy and Lyα halo) to be identical by considering normal distribution. We consider a Lyα halo as detected if p₀ ≤ 0.05. Figure 5 shows p₀ as a function of Lyα luminosity. Out of the 184 objects for which we have a Lyα halo scale length measurement, 20 do not show statistical evidence for extended Lyα emission, mainly due to the large errors on their size measurements. Interestingly, some of these are relatively bright in Lyα (e.g. #218 in Fig. 3). Thirty objects have a Lyα halo with a very high significance (p₀ ≤ 10^-5).

Finally, out of the objects for which we have reliable Lyα halo measurements and excluding the objects with only upper limits to their Lyα halo scale length, 145 galaxies have a statistically significant Lyα halo (116 and 29 from the mosaic and udf-10 field, respectively), which represents ≈80% of the sample (145 out of 184 objects).

5.2. Halo sizes

Figure 6 shows the distribution of Lyα halo scale lengths for the 145 haloes that are considered to be statistically significant. This distribution contains scale lengths ranging from 1.0 to 18.7 kpc with a median value of 4.5 kpc. For reference, Xue et al. (2017; hereafter X17) measured similar halo scale lengths as our median scale length value; i.e. 5–6 kpc from their median stacking images of all LAEs. The W16 authors measured slightly smaller scale lengths, with a median value of 3.4 kpc, but still in good agreement with our results.

Figure 6 also shows an extended tail of large halo scale lengths (>7 kpc). These large haloes represent less than 12% (29 galaxies) of our total sample, plausibly explaining why they are not present in the smaller W16 sample.

5.3. Fluxes and equivalent widths

5.3.1. Lyα flux and equivalent width fractions

Our 2D, two-component decomposition of the Lyα spatial distribution provides estimates of the decomposed flux from the core, F_cont, and from the halo, F_halo. The median halo flux is ≈5.4 × 10^-18$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$ and ≈9.4 × 10^-18$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$ for the udf-10 and mosaic data cubes, respectively. Following the same approach as W16, we defined the Lyα halo flux fraction, X_Lyα,halo, as F_halo/ (F_halo + F_cont), quantifying the contribution of the halo to the total Lyα flux. For galaxies with a UV continuum detection and with reliable Lyα halo scale length measurements, we find halo flux fractions ranging from ≈17 to ≈99% (with an average of 65%). Our results are in very good agreement with W16 (they find a mean value of 70%), although we find more haloes with small Lyα halo flux fractions than W16. This difference is likely due to our larger sample and to the improved quality of the UDF data cubes with respect to the HDFS (see Figs. 8 of B17).

Figure 7 shows the Lyα halo flux fraction as a function of the core and halo scale lengths and as a function of the total Lyα luminosity measured by the CoG method (see Sect. 5.3.2). In order to test for correlations, we calculated the Spearman rank correlation coefficients ρ_s (e.g. Wall 1996) and the corresponding p-values p₀, which correspond to the probability of the null hypothesis that no monotonic relation exists between the two variables. The value ρ_s varies between − 1 and + 1, with 0 implying no correlation. The test does not take the error bars into account. We obtained (ρ_s = −0.05, p₀ = 0.52), (ρ_s = 0.03, p₀ = 0.69), and (ρ_s = −0.15, p₀ = 0.05) for the X_Lyα,halo−rs_halo, X_Lyα,halo−rs_cont, and X_Lyα,halo−L_Lyα, respectively. We find no correlation between the fraction of Lyα emission in the halo and (i) the Lyα halo scale lengths (second panel); (ii) the UV continuum scale lengths (third panel); and (iii) the total Lyα luminosities (right panel).

Figure 8 shows the Lyα flux in the halo as a function of its halo scale length. Our limiting Lyα halo flux and scale length (see Sect. 4.3.2) are indicated. We find a clear correlation between these two properties (ρ_s = 0.4, p₀< 10^-8). While this correlation is partially created by our halo flux limit (red lines), the correlation is still readily apparent for brighter haloes. The correlation is positive, such that Lyα haloes with larger scale lengths have higher halo fluxes.

5.3.2. Total Lyα flux and equivalent width

The Lyα flux was computed by integrating inside the circular aperture corresponding to the CoG radius. This radius (r_CoG) was determined by averaging the flux in successive annuli of 1 pixel thickness around the Lyα emission centre until a certain annulus for which the averaged flux reaches the noise value. The centre of this last annulus corresponds to r_CoG. From this aperture, we extracted a spectrum and integrated the flux corresponding to the Lyα line width; the borders of the line are set when the flux goes under zero. These spectral bandwidths are indicated by the vertical black dotted lines in the third panel of Figs. 2 and 3¹.

This method ensures that most of the Lyα flux is encompassed for each object, which is not the case if we use a single fixed aperture for all of the objects. Figure 1 (lower panel) shows the distribution of total Lyα fluxes for our sample, which spans 2 orders of magnitude from 1.74 × 10^-18$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$ to 1.12 × 10^-16$\mathrm{ergs}{-1}^{}\mathrm{cm}{-2}^{}$.

Rest-frame Lyα equivalent widths (EWs) were calculated using the UV continuum measured by Hashimoto et al. (2017; hereafter H17). The H17 authors performed careful UV continuum measurements using several HST bands (2 or 3 bands depending on the object). After cross-matching the respective catalogues, we obtained Lyα EW measurements for 155 of the 184 galaxies for which we have a halo scale length measurement.

Figure 9 shows the distribution of Lyα halo scale lengths as a function of total Lyα luminosity and rest-frame Lyα EW. With respect to previous studies that employed stacking, our sample goes much deeper and we probe much smaller Lyα haloes. The Spearman rank correlation test provides ρ_s = 0.222 (p₀ = 0.002) for the rs_halo−L_Lyα relation and ρ_s = 0.09 (p₀ = 0.29) for the rs_halo−EW₀ relation.

We find a suggestion of a correlation between the Lyα scale length and total Lyα luminosity, albeit with very large scatter. In particular the bright LAEs tend to have large haloes (rs_halo ≳ 3 kpc), whereas there is more dispersion at lower Lyα luminosities. The X17 authors found a clear correlation in their stacks corresponding to our bright LAEs whereas W16 found no such correlation. By computing the Spearman rank correlation coefficients for our objects in their luminosity range (41.6 < log (L_Lyα) < 42.7), we find no correlation either (ρ_s = 0.07, p₀ = 0.59).

We also do not find a correlation between halo sizes and rest-frame EWs. Both W16 and X17 found a similar result. As an aside, it is worth noting that some objects have very large EWs (exceeding 200 Å). The H17 work provides for a detailed analysis of these objects.

5.4. Lyα line profiles

Next, we explore the connection between the spectral and spatial properties of the Lyα emission. The diversity of Lyα line profiles can be appreciated by looking at Figs. 2 and 3. While most of the lines are asymmetric and single peaked, others appear to be double peaked (see objects #1087 and #106 of Figs. 2 and 3, respectively). This diversity is directly reflected in the Lyα FWHM measurements, which span a large range from 118 to 512 km s^-1.

We performed the measurement of the Lyα FWHM only on the single-peaked Lyα lines so that objects with doubled-peaked profiles are excluded. If the blue peak is comparable in flux to the red peak, the Lyα line is referred to as a Lyα doublet (see object #106 in Fig. 3), whereas if the blue peak is much fainter, the feature is referred as a blue bump (see object #1087 in Fig. 2). We carried out an inventory of the various line profiles encountered in our sample. Out of our 252 galaxies, 15 objects show a Lyα line with a blue bump and 8 objects have a Lyα doublet. Put together, the double-peaked profiles therefore represent a small fraction (<10%) of our sample. The halo properties of such double-peaked line objects are not significantly different from those of the rest of the sample.

We connect the FWHM of single-peaked Lyα lines with Lyα halo sizes in Fig. 10. The smallest Lyα haloes, for which we only have an upper limit on their scale length or halo flux (see Sect. 4.3.2) and the Lyα haloes with low statistical significance (see Sect. 5.1) appear to have a narrower Lyα line (<350 km s^-1), whereas the galaxies with significant extended Lyα emission span a wider range of FWHM values.

Figure 11 shows the Lyα EW plotted against the FWHM, colour coded by the Lyα halo flux fraction. We show in this figure only the 121 objects with a statistically significant Lyα halo (see Sect. 5.1) and an EW measurement (see Sect. 5.3.2). It is also apparent that we do not find evidence for a significant anti-correlation between the EWs and the FWHMs of the Lyα lines (ρ_s = −0.21, p₀ = 0.02).

The objects that have less than 30% of their total Lyα flux in the halo appear to have narrower Lyα lines than the rest of the sample (<300 km s^-1). Certainly, the objects with a large Lyα line width (>400 km s^-1) have small EW (<100 Å) and >50% of the Lyα flux is in the halo.

6. Connecting host galaxies to Lyα haloes

In this section, we investigate the connection between Lyα properties and the properties of the host galaxies. First we consider the general UV properties. We then compare galaxy and Lyα halo sizes. Finally, we explore the coevolution of UV and Lyα halo sizes with redshift.

6.1. UV properties

The upper panel of Fig. 12 shows the expected correlation (ρ_s = 0.35, p₀ ~ 10^-6) between the UV sizes and UV magnitudes of galaxies (Shibuya et al. 2015). The lower panel shows Lyα halo scale length as a function of absolute far-UV magnitude. According to the Spearman test coefficient (ρ_s = 0.19, p₀ = 0.02), there is a suggestion of a positive correlation between the Lyα halo size and UV magnitude for our selected sample of LAEs albeit the scatter is large. If the correlation is real, it would agree with the results of X17 who found that Lyα halo sizes are positively correlated with UV luminosities. Similarly, W16 found that UV-luminous galaxies (M_UV< −19) tend to have Lyα haloes with larger scale lengths (rs_halo ≳ 3 kpc). This result is also observed in our larger sample.

6.2. Comparison of sizes

Next we compare the UV continuum and Lyα emission scale lengths resulting from our 2D two-component model. Figure 13 shows the Lyα halo scale lengths plotted as a function of UV continuum scale length. First, according to the Spearman correlation test, Lyα scale lengths are positively correlated (ρ_s = 0.32, p₀ ~ 10^-5) with galaxy UV sizes (albeit with large scatter). Indeed, the Lyα scale lengths are always between ≈4 and >20 times larger than the continuum scale length with a median size ratio of 10.8 (the lower quartile at 6.0 and the upper percentile at 19.1). This results is in very good agreement with W16, albeit this scatter in the ratio of scale lengths exceeds the range from their sample. This plot also shows that we do not detect any LAEs without a Lyα halo. This is valid for 145 galaxies (80%) of our sample given that the remaining 39 galaxies (20%) only have either upper limits (grey area) or large error bars on their halo scale lengths (empty circles in Fig. 13 upper panel). As such, we can only confirm the absence of a Lyα halo around a galaxy larger than this size limit (i.e. rs_cont ≳ 1 kpc; see Sect. 4.3.2). This condition is shown as the green area in Fig. 13. The dashed areas show the range of our wavelength-dependent detection size limit (see Sect. 4.3.2).

6.3. Size evolution

The evolution of both the UV and Lyα halo sizes is shown in Fig. 14. While the UV size of galaxies decreases with redshift as expected (Shibuya et al. 2015), we do not find that Lyα halo sizes show significant evolution between redshifts 3 and 6. The W16 authors found that Lyα halo sizes decrease with increasing redshift. However, their sample at z> 5 consists only of five galaxies and the dispersion is large. Because we do not have enough objects in the higher redshift bin, we conclude that there is currently no clear evidence supporting an evolution of Lyα halo sizes with redshift above z = 3. Momose et al. (2014) also investigated the size evolution of their stacked LAEs, finding no evidence for evolution of Lyα halo sizes from z = 2.2 to z = 5.7 and a possible but very uncertain increase from z = 5.7 to z = 6.6. This result implies a higher Lyα/UV scale length ratio at high redshift and hence suggests that the fraction of CGM probed by the Lyα emission is actually increasing with redshift as the galaxies are known to be more compact and less massive as high redshift (Shibuya et al. 2015).

We now compare our scale length measurements with 12 local starburst galaxies (0.028 <z< 0.18) from the LARS sample (Hayes et al. 2013; Guaita et al. 2015). These authors measured the spatial extent of the Lyα emission using the Petrosian 20 percent radius R_p20 (Petrosian 1976) and compared to the corresponding radius measured from Hα. Similar to the high-redshift galaxies, some local galaxies show extended Lyα emission (seven galaxies according Hayes et al. 2014). Their resulting Lyα/Hα size ratios range from 1 to 3.6, with an average of 2. For comparison, we also calculate the Petrosian radii of our galaxies. As in W16 (see their Fig. 12), our galaxies at z> 3 appear to have Lyα haloes with larger Petrosian radii as well as higher Lyα/UV size ratios than local galaxies.

7. Discussion

7.1. Probing the CGM

7.1.1. Ubiquity of Lyα haloes around LAEs

The high fraction of LAEs with a significant Lyα halo (≈80%) demonstrates that Lyα haloes are a common property of LAEs at high redshift. As such, this also suggests that the CGM has a rich “cold” gas content. Theoretical analyses and numerical simulations indeed predict that neutral hydrogen should be present around high-redshift star-forming galaxies (e.g. Kereš et al. 2005; Fumagalli et al. 2011; Faucher-Giguère et al. 2015). Our results therefore appear consistent with the canonical vision of the galaxy formation at high-z. For the remaining 20% we have only upper limits or very uncertain Lyα halo size measurements, which prevents us from drawing firm conclusions.

That such a result is not seen for local galaxies (cf. LARS sample; Östlin et al. 2014; Hayes et al. 2013, see Sect. 6.3) underlines that there is a clear evolution of the CGM across cosmic time. In the LARS sample, Lyα haloes are only detected in 50% of cases (Hayes et al. 2014). Moreover, Hayes et al. (2013) found the Lyα emission to be more extended than the UV continuum by a factor of 2.4 on average (interstellar medium scales), whereas we find a factor of ≈10 for our sample (CGM scales). This difference could be due to an evolution of the Lyα escape fraction with redshift, possibly due to dust content evolution (Hayes et al. 2011; Dijkstra & Jeeson-Daniel 2013) or because of sensitivity limitations. This can also suggest that the contribution of the mechanisms powering the Lyα haloes evolves with cosmic time.

This study is limited to Lyα emitters. In future work, we intend to search for extended Lyα emission around individual galaxies that were not detected based on their Lyα line (i.e. LBGs). For example, Steidel et al. (2011) detected extended Lyα emission around a stacked Lyα absorber galaxies sample and around massive LBGs. This promising result motivates us to search for Lyα haloes around all high-redshift galaxies with MUSE.

7.1.2. Lyα halo size – host galaxy correlations

Supposing Lyα haloes are a general property of star-forming galaxies, it is interesting to know how the halo sizes correlate with other properties of the host galaxies. We searched for such correlations in Sect. 6. Both the Lyα halo flux and scale length seem to correlate with both UV magnitude and scale length, suggesting that the Lyα halo properties are actually linked to the UV properties of the host galaxy. Consequently, this may suggest that the star formation rate directly influences the powering of the Lyα haloes. Also, if the trends are real, the correlations suggest that the UV-brighter objects are associated with different physical conditions, such as kinematics, gas content, and distribution, which favour the production of Lyα haloes compared to the fainter objects from our sample.

7.1.3. Lyα spatial extent versus virial radius of DM haloes

We now attempt to assess the maximum CGM scales that are traced by Lyα emission as a function of UV magnitude. To do so, we compare the maximum detected extent of the Lyα haloes (measured using the CoG method) with the virial radius of the dark matter (DM) haloes of galaxies predicted by the semi-analytic model of Garel et al. (2015; see contours in the upper panels of Fig. 15). We compare those two extents for four different redshift bins (z ≃ [3,4,5,6]). In the lower panel, the purple dotted lines indicate the median value of the ratio of CoG radii over the mean virial radii r_CoG/ ⟨r_vir⟩ in each redshift bin. This ratio appears to increase with redshift (ratio median values of [57%, 64%, 69%, 87%] for median redshift bins of [3.2, 3.8, 4.8, 5.9]) suggesting that Lyα emission is probing a larger percentage of the CGM towards high redshift. This is not surprising because our measured Lyα halo sizes do not show any evolution with redshift while the Garel et al. (2015) model predicts that, at fixed M_UV, galaxies reside in smaller DM haloes at higher redshift.

While a weak anti-correlation can be guessed, the fraction of CGM probed by the Lyα emission is fairly constant with UV magnitude in each redshift bin.

Our Lyα haloes therefore reach on average more than ≳50% of the predicted virial radius of their host galaxy (irrespective of M_UV) and go even beyond for higher redshift. This result demonstrates that the Lyα emission is a powerful tracer of the gas located inside the virial radius (e.g. the CGM) but not at larger scales (e.g. IGM) considering our current detection capacities.

7.2. Origin of the Lyα haloes

Taking advantage of our large statistical sample, we now attempt to assess the contribution of the different proposed Lyα emission processes that could be responsible for our observed Lyα haloes. For each process, we review the emission mechanism, discuss the expected observational signatures and, where possible, compare these expectations with our results.

7.2.1. Stellar origin with scattering in an outflowing medium

The scattering of Lyα photons produced in star-forming regions is one of the candidates to explain Lyα haloes. For this mechanism, Lyα photons are produced by recombination associated with the stellar UV radiation in the HII regions of galaxies. A fraction of those Lyα photons can be absorbed by interstellar dust but the escaping photons scatter into the surrounding neutral hydrogen gas and can be redirected towards the observer, leading to the observed Lyα haloes.

The main question is therefore whether the stellar content of the galaxies produces enough ionizing photons to power the observed Lyα emission. In a similar approach to W16 (see their Sect. 7.2 for more details), the condition to be tested is a condition on the Lyα EW as this quantity gives a direct comparison between the continuum and Lyα fluxes. The maximum dust-free Lyα EW estimated for a stellar origin ranges from ≈50 to 200 Å (Charlot & Fall 1993). While ≈17% of our sample has Lyα EWs higher than 200 Å, which suggests that Lyα photons do not only come from the HII regions, most of our galaxies have Lyα EWs lower than 200 Å. This suggests that the stellar UV continuum alone can power the haloes. In any case, given that the EW depends on stellar metallicity and initial mass function and can be affected by bursty star formation histories (Schaerer 2003; Raiter et al. 2010), the objects with EW> 200Å values may be interpreted without invoking other Lyα production channels. A more detailed discussion of the objects in our sample with large Lyα EWs is presented in H17.

Information is also encoded in the spectral shape of the Lyα line. Looking at our sample, most of our Lyα spectra show a single asymmetric line, as expected for Lyα scattering processes in outflowing media (Verhamme et al. 2006; Dijkstra & Kramer 2012; Yang et al. 2016). Outflows facilitate the escape of Lyα photons emitted in star-forming regions from the ISM (Dijkstra et al. 2006; Verhamme et al. 2012; Behrens et al. 2014; Behrens & Braun 2014) and can therefore be responsible for the observed Lyα haloes. Some observational evidence has been found supporting this scenario in local galaxies (e.g. Bik et al. 2015; Herenz et al. 2016). However, the Lyα spectra do not indicate where the Lyα photons are produced. Indeed, they can be produced either in HII regions well within the galaxy and then scatter in the CGM or ISM (some very compact LAEs do show asymmetric Lyα line profiles) or in the CGM and still scatter within the CGM producing asymmetric lines (Cantalupo et al. 2005). Analyses of spatially resolved spectra along with Lyα transfer simulations should however be able to help disentangle between the different effects.

If the observed Lyα haloes are powered by Lyα radiation produced inside the galaxies and scattered outwards, we expect the spatial and spectral properties of these haloes to correlate (Verhamme et al., in prep.). In particular, the halo flux fraction is predicted to increase with the spectral shift of the peak and the FWHM of the Lyα line. Figure 11 however does not show such a trend.

Hence, our results indicate that the scattering of Lyα photons created in HII regions can contribute to the powering of Lyα haloes but it is difficult to quantify their contribution.

7.2.2. Gravitational cooling radiation

Another scenario to explain the extended Lyα emission is the so-called “cooling radiation” (e.g. Haiman et al. 2000; Fardal et al. 2001; Furlanetto et al. 2005). In this process, Lyα photons are emitted by collisionally excited circum-galactic gas, which converts gravitational energy into kinetic and thermal energy as it falls into the DM halo potential. Cooling radiation has been postulated to come into play at large radii, where the Lyα photons are less likely to be absorbed by dust (Fynbo et al. 2001; Fardal et al. 2001). However, because the density is higher at the centre of the DM halo and the Lyα emissivity resulting from cooling is proportional to the density squared if the gas is warm enough, the Lyα emission is expected to be centrally concentrated (as found in Rosdahl & Blaizot 2012). The bulk of the Lyα radiation that we would observe from such a geometry would therefore have scattered outwards through an infalling scattering medium.

A number of other theoretical analyses have been carried out to predict the expected cooling radiation contribution (Furlanetto et al. 2005; Dijkstra & Loeb 2009; Faucher-Giguère et al. 2010; Rosdahl & Blaizot 2012). Such numerical simulations are nevertheless difficult to perform as they require high resolution and expensive radiative transfer treatments. Recently, Lake et al. (2015) performed hydrodynamic and radiative transfer simulations of LAEs and found that star formation accounts for the origin of the majority of diffuse Lyα emission but that cooling radiation can also have a significant contribution; i.e. 40–55% of the total Lyα luminosity within distances up to the LAE virial radius.

According to theoretical predictions (Dijkstra et al. 2006) and confirmed by numerical experiments (Trebitsch et al. 2016), for Lyα radiative transfer in a cooling gas, the resulting Lyα line is expected to be blueshifted with a blue tail with respect to the centre of the line. The effect of IGM absorption, however, even at the redshifts considered here, can have a significant impact on the blue side of the Lyα line profile (Laursen et al. 2011). Most of our objects do not show any of the non-resonant lines needed to determine the systemic redshift of the galaxy and so to precisely measure such a blue shift. However, looking at our sample, none of our LAEs shows a clear single-peaked asymmetric line towards the blue. The predicted blueshifted feature could be manifested as the blue bumps observed for 10% of our LAEs. Such lines indeed show a shifted blue peak and an enhanced red peak, suggesting that cooling radiation cannot fully account for the shape of the Lyα lines. However, some theoretical predictions are obtained by averaging over all directions, which is of course not the case for observed spectra. Moreover, the line of sight can strongly impact the line profiles (Verhamme et al. 2012; Gronke & Dijkstra 2014). As such, the spectral features of the Lyα line should be interpreted with caution here.

Alongside Lyα spectral properties, Lyα luminosities also provide crucial information. Both Rosdahl & Blaizot (2012) and Dijkstra & Loeb (2009) predict the Lyα luminosity produced by cooling radiation in a 10¹¹M_⊙ DM halo to be 5 × 10⁴¹ erg s^-1. As our Lyα halo luminosities (L_Lya,halo) are higher, this suggests that cooling is not the only process producing the Lyα halo emission or that the bulk of our LAEs reside in DM haloes more massive than ≈ 10¹¹M_⊙. This latter option is unlikely, however. In Fig. 16, we plot the UV luminosity (L_UV) to Lyα halo luminosity ratio as a function of L_UV for the sources in our sample². First, the anti-correlation between L_Lyα,halo/L_UV and L_UV shows that the halo component contributes more in UV-faint galaxies than in brighter UV sources. Interestingly, this trend may actually reflect the so-called “Ando effect”, i.e. the fact that faint M_UV objects appear to have large Lyα EW³, which is commonly observed at high redshift (Ando et al. 2006, H17). Second, we perform a robust linear fit to the data and we measure a slope of − 0.52 ± 0.05 for the L_Lyα,halo/L_UV versus L_UV relation. This corresponds to ${\mathit{L}}_{\mathrm{Ly}\mathit{\alpha ,}\mathrm{halo}}\mathrm{\propto }{\mathit{L}}_{\mathrm{UV}}^{\mathrm{0.45}}$ and denotes that UV bright galaxies in our sample have more luminous L_Lyα,halo haloes. This result can be directly compared to the predictions of Rosdahl & Blaizot (2012) who only considered cooling radiation at z = 3. These authors found a slope of 0.625 (if we assume UV luminosity proportional to the square of the DM halo mass), which is different from our result but not so dissimilar. Hence we cannot rule out this scenario. In the “scattering from HII regions” scenario, one would expect that the flux in the Lyα halo would globally scale with the number of Lyα photons that escape the galaxy, i.e. L_UV times the Lyα escape fraction from the ISM. We should therefore observe L_Lya,halo ∝ L_UV; if all galaxies have the same Lyα escape fraction and if there are more neutral hydrogen atoms in the CGM than Lyα photons from the galaxies. It is worth noting however that varying dust content or ISM column density could strongly affect this relation by introducing a UV magnitude-dependent Lyα escape fraction.

Put together, our analysis does not allow us to give a firm conclusion about the contribution of cooling radiation in the production of Lyα haloes.

7.2.3. Lyα fluorescence

Another possible origin for the Lyα haloes is the Lyα fluorescence resulting from the recombinations of hydrogen that is photo-ionized by Lyman continuum (LyC) radiation generated by nearby quasars, young stars, or by the cosmic UV background (UVB) (Furlanetto et al. 2005; Cantalupo et al. 2005; Kollmeier et al. 2010). This scenario is usually invoked for giant Lyα nebulae, within which quasars are known to reside (Cantalupo et al. 2014; Borisova et al. 2016; and see Cantalupo 2017, for a review), as well as for compact dark galaxy sources (Cantalupo et al. 2012; Fumagalli et al. 2016; Marino et al. 2017).

According to the predictions of Haardt & Madau (1996) and Cantalupo et al. (2005), the resulting Lyα SB produced by the diffuse ionizing background is significantly lower (~10^-20 erg s^-1 cm^-2 arcsec^-2 at z ≈ 3). The expected effects of the UVB therefore appear to be negligible for the individual objects of our study.

According to the calculations of Gallego et al. (2017), which uses the same MUSE UDF data as our study, the required LyC escape fraction from the ISM for stars to produce the observed Lyα halo in their stack of LAE pairs (SB of ~ 3 × 10^-20 erg s^-1 cm^-2 arcsec^-2) is extremely small (f_esc ≈ 0.02). This result suggests that it is not so difficult to have a large ionized fraction of gas in the inner parts of the haloes; the high gas densities and clumpiness of the medium moreover favour the Lyα fluorescence.

Consequently, we cannot rule out the contribution of the Lyα fluorescence for the powering of our observed Lyα haloes.

7.2.4. Satellite galaxies

Momose et al. (2016) proposed another explanation where the Lyα halo would be powered by several satellite galaxies emitting Lyα emission around the central galaxy. Shimizu & Umemura (2010) and Lake et al. (2015) have shown using cosmological simulations that Lyα haloes are indeed associated with such surrounding galaxies. Recently, Mas-Ribas et al. (2017) applied an analytic formalism (Mas-Ribas & Dijkstra 2016) to investigate the plausibility of this scenario by using various satellite clustering conditions. These authors found that satellite sources can indeed play a role in the powering of Lyα haloes at large distances (20 ≲ r ≲ 40 physical kpc) from the galaxies. According to their modelling, such satellite galaxies would be very faint in the UV continuum (M_UV> −17) so that they would be undetectable by any current instruments and may therefore be missed in current surveys.

Applied to the case of our data, we can expect that in the presence of satellites, which emit Lyα emission and are undetected in UV in the HST images, our Lyα haloes would appear clumpy and rather asymmetric. We do not observe such clumpiness in the central regions of our Lyα NB images. However the MUSE PSF acts to smooth out Lyα clumps, making their detection impossible. Furthermore, at larger radii the S/N of the Lyα NB image drops significantly, making the detection of clumps or asymmetries very challenging. Moreover, a large contribution from the star formation in satellites is expected to provide similar UV and Lyα extended emission. Our measured UV continuum scale lengths however appear much smaller than the Lyα scale lengths (see Sect. 6.2). In the presence of Lyα-emitting satellites we might also expect the Lyα haloes would be offset from the host galaxy. Such an offset is not observed for most of our objects.

Put together with all these elements, it seems somewhat unlikely that there is a significant contribution from satellite galaxies to the powering of Lyα haloes. We cannot however completely rule out the possibility that unidentified satellites partly power the Lyα haloes.

7.2.5. Future directions

In summary, our results are suggestive of a scenario, in which the following range of processes can be responsible for the observed Lyα haloes:

• The scattering on CGM scales of Lyα photons that are produced by recombinations in HII regions.

• The cooling radiation triggered by gas inflowing onto the host galaxies.

• The Lyα fluorescence associated with hydrogen recombinations after ionization by Lyman continuum (LyC) radiation present in the CGM.

While those processes have to be considered together, their respective contributions cannot be constrained by our data. To try to disentangle the relative impact of the different processes, we need to know where the Lyα photons are produced, which is not straightforward because Lyα is a resonant line. As such, more observations are needed. In particular, Hα observations by the James Web Space Telescope (JWST) will directly tell us the origin of the Lyα emission as it is not a resonant line. Hα emission that is more extended than the UV continuum would be a direct piece of evidence that the Lyα emission is produced in the CGM (i.e. by fluorescence). On the other hand, compact Hα emission would indicate that Lyα photons are produced in the ISM and then propagate in the CGM by resonant scattering. The adaptive optics (AO) on MUSE, currently in commissioning, also promises good progress as it will significantly improve the PSF and therefore enable the detection of smaller haloes and allow a precise characterization of the halo morphologies. Finally, the help of theoretical and numerical studies will be needed to fully understand the processes at play and their respective contributions.

8. Summary and conclusions

Thanks to the significant increase in sensitivity enabled by the MUSE instrument, we studied the CGM gas content of an unprecedentedly large sample of individual star-forming galaxies at redshift z = [ 3–6] in the Hubble Ultra Deep Field. Our LAE sample was selected to have a good S/N (>6) and to be isolated (see Sect. 2.2). Our galaxy-by-galaxy based analysis allows us to characterize individual Lyα halo properties and to explore possible correlations with the UV properties of the host galaxies. Our major results are summarized as follows:

• 1.

  We detect diffuse Lyα emission with high confidence around 145 individual LAEs. This represents 80% of our objects for which we have reliable Lyα halo measurements. Among the objects for which we have reliable Lyα halo scale length measurements, 20 do not show a significant Lyα halo detection, mainly owing to large errors on their halo size measurement (see Sect. 5.1). Put together, our data suggest that extended Lyα haloes are common around Lyα emitters at high redshift.

• 2.

  We find a large range of Lyα halo scale lengths, emphasizing the diversity of configurations of the cool CGM. The halo scale lengths in our sample range from 1.0 to 18.7 kpc with a median value of ≈4–5 kpc. We also show that the Lyα emission probes the CGM out to large radii (Fig. 15), reaching on average ≈50% of the virial radius according to the comparison of our data with predictions from a semi-analytic model. This result shows that Lyα emission is a powerful tool to map the cold hydrogen around high-redshift galaxies.

• 3.

  The Lyα haloes properties of our selected sample of LAEs appear to be dependent on the stellar content of the galaxies. Both Lyα halo spatial extents and fluxes are found to be positively correlated with UV magnitudes and spatial extents of the host galaxies, although the correlation with UV magnitude is not as clear.

• 4.

  While Lyα halo scale lengths appear to be considerably larger at z> 3 than at z ≃ 0 (from a comparison with the LARS sample), we do not observe any significant evolution of the Lyα scale lengths between redshift 3 and 6. This implies an evolution of the CGM content between z ≃ 0 and z = 3.

• 5.

  The galaxies that are less likely to have a Lyα halo as well as those with small haloes (rs_halo ≲ 1 kpc, i.e. objects with upper limits on the halo size), have on average narrower Lyα lines than the rest of the sample (Fig 10). This suggests that the Lyα line is less broadened when the gas content in the CGM is low. However, it is worth noting that Lyα line of galaxies with a high confidence Lyα halo are not systematically broader.

• 6.

  With the information from our data we attempt to explore the origin of the Lyα haloes around star-forming galaxies. While we find no evidence for a dominant contribution from a single particular process, we are not able to rule out any of the scenarios we consider, i.e. scattering from star-forming regions, fluorescence, and cooling radiation from cold gas accretion, except maybe the scenario for which satellites significantly contribute. Indeed, while we do not find significant evidence for the “satellite scenario”, our data cannot disentangle whether the Lyα photons are produced in the star-forming regions and then scatter in the CGM or “in-situ” in the CGM from gravitational cooling radiation and/or from fluorescence. As a consequence, further observations and analysis will be needed to understand the powering process(es) of the Lyα haloes (JWST, MUSE with AO, and theoretical and numerical analyses).

The MUSE instrument has enabled us to extend the sample of measurements of Lyα haloes to fainter and smaller galaxies, which are more representative of the bulk of the galaxy population. Our study underlines the significant cold gas content of the Universe between redshifts 3 and 6, regardless of the nature of the Lyα halo emission mechanism.

This new study highlights that Lyα emission presents an exciting new opportunity to study the diffuse and low-SB gas in the vicinity of faint high-redshift galaxies. In the coming years, adaptive optics mounted on MUSE/VLT will allow us to be even more precise in the detection and characterization of these Lyα haloes. By improving the PSF, it will be possible to detect smaller Lyα haloes and thus confirm if there are LAEs without any halo component. Within the context of the substantial recent progress in improving numerical simulations and theoretical models, detailed comparisons of models and observations of Lyα haloes around normal star-forming galaxies are now possible and promise to significantly expand our understanding of the mechanisms that regulate the gas that flows in and out of galaxies in the early Universe.

---

Acknowledgments

F.L., R.B., S.C., H.I., and M.A. acknowledge support from the ERC advanced grant 339659-MUSICOS. TG is grateful to the LABEX Lyon Institute of Origins (ANR-10-LABX-0066) of the Université de Lyon for its financial support within the programme “Investissements d’Avenir” (ANR-11-IDEX-0007) of the French government operated by the National Research Agency (ANR). R.B., T.C., and J.R. acknowledge support from the FOGHAR Project with ANR Grant ANR-13-BS05-0010. S.C. acknowledges support from Swiss National Science Foundation grant PP00P2 163824. T.C. acknowledges support from the OCEVU Labex (ANR-11-LABX-0060) and the A*MIDEX project (ANR-11-IDEX-0001-02) funded by the “Investissements d’avenir” French government programme managed by the ANR. J.R. acknowledges support from the ERC starting grant CALENDS. J.S. acknowledges support from ERC Grant agreement 278594-GasAroundGalaxies. R.A.M. acknowledges support by the Swiss National Science Foundation. J.B. acknowledges support by Fundação para a Ciência e a Tecnologia (FCT) through national funds (UID/FIS/04434/2013) and Investigador FCT contract IF/01654/2014/CP1215/CT0003, and by FEDER through COMPETE2020 (POCI-01-0145-FEDER-007672).

References

Appendix A: Effect of the exposure time on the Lyα halo detection

Our sample is constructed from two overlapping data sets with various exposure times, where the udf-10 data is 3 times deeper on average than the mosaic data. We are therefore able to investigate how our Lyα halo size measurements vary with depth.

Within our sample, 26 objects are both detected in the udf-10 and mosaic fields but only 15 have a reliable Lyα halo measurement (see Sect. 4.3.2). Figure A.1 shows the difference of the Lyα halo scale length measurements against the S/N (left panel) and the total Lyα flux (right panel) of the Lyα NB image constructed from the mosaic data cube for these 15 objects. The median difference is small (<0.1 kpc) and the error bars encompass the measured offsets for every object.

Figure A.2 shows a more detailed comparison of the two data cubes for three representative objects. The S/N of Lyα NB images from the udf-10 is on average larger by a factor ≈2 compared to the mosaic. This is consistent with the noise analysis of the UDF MUSE data cubes given in B17. This is also clearly highlighted by the contours, which are clearer and less splintered in the udf-10 images. The NB images are optimized in terms of S/N from both data cubes. Hence, the spectral bandwidths are not always the same in the two data cubes for a given object. As a consequence the NB images and therefore the Lyα centroid measurement can be slightly different explaining that the SB profiles do not perfectly overlap in the inner region.

The first object (#180, top row) shows a discrepancy in the halo scale length measurements between the two data cubes (object indicated as red point in Fig. A.1). Visual inspection of radial SB profiles (right panel) shows that the Lyα halo of this object is lost in the noise for the mosaic data cube. We apply a S/N cut of 6 to define our sample (see Sect. 4.3.1) The S/N of the object we are showing here is higher (9.3) than the S/N cut (6) and this poses questions about the reliability of our S/N cut.

The next object #168 (middle row of Fig. A.2, orange dot in Fig. A.1) offers a counter example. While the S/N of the mosaic Lyα NB image is lower than the previous example and very close to our S/N cut value (6.1), the scale lengths measured in the two data cubes are similar for this example. These two examples highlight that the S/N cut we define using simulated extended objects is not absolute and that for S/N < 10, halo sizes can be underestimated in some cases. This uncertainty is however encompassed in the error bar.

The last example (#149, last row of Fig. A.2 and green dots in Fig. A.1) shows the comparison of the detection for an object with a good S/N in the two data cubes. Reassuringly, the scale lengths measured in the NB images are similar with a larger error on the scale length fit to the mosaic data. Put together, this illustrates the importance of surface brightness sensitivity for the detection of extended Lyα emission.

Appendix B: Table of Lyα halo measurements

All Tables

All Figures