Kinematics, turbulence, and star formation of z ∼ 1 strongly lensed galaxies seen with MUSE Free

Abstract

1 INTRODUCTION

High-redshift disc galaxies display some striking differences when compared with their local counterparts: not only do they harbour giant H ii star-forming regions, but they also have higher gas velocity dispersions and higher gas fractions than local discs (see Glazebrook 2013 for a review). The star-forming regions seen in these high-z discs, referred to as clumps or knots (usually identified in rest-frame UV/optical images; e.g. Elmegreen & Elmegreen 2006), are also more extreme than the local star-forming regions. They can have sizes of up to one kiloparsec, star formation rates (SFRs) between 0.5 and 100 M_⊙ yr⁻¹ and masses up to 10⁸–10¹⁰ M_⊙ (e.g. Swinbank et al. 2009; Jones et al. 2010; Dessauges-Zavadsky et al. 2011; Förster Schreiber et al. 2011; Guo et al. 2012), which makes them significantly more massive than giant clouds in local galaxies, with masses of 10⁵–10⁶ M_⊙. While it was initially speculated that these features had their origin in merging episodes, the advent of Integral Field Spectrographs showed that most clumpy galaxies have smooth velocity fields and are rotationally supported (e.g. Genzel et al. 2006; Bouché et al. 2007; Cresci et al. 2009; Förster Schreiber et al. 2009; Wisnioski et al. 2015; Contini et al. 2016), suggesting that these clumps may be part of the secular evolution of disc galaxies.

The physical picture that explains why rotating discs are more clumpy and turbulent at z ≥ 1 is still being investigated, both observationally and theoretically. One possible scenario is that the highly turbulent interstellar medium (ISM) is stirred by radiation pressure and winds of strong star formation taking place in these young galaxies (e.g. Lehnert et al. 2009; Green et al. 2010; Lehnert et al. 2013). Another possibility is that the inflow of gas into the galaxy from the cosmic web provides enough energy to sustain these high-velocity dispersions, although some works seem to point that this energy is insufficient to maintain the high-velocity dispersions over long time-scales (e.g. Elmegreen & Burkert 2010). Finally, a third scenario is gravitational instability: the high gas fractions of these galaxies make them gravitationally unstable, causing gas to spiral down to the centre of the galaxies, converting gravitational energy into turbulent motions and driving galaxies to stability (e.g. Bournaud, Elmegreen & Elmegreen 2007; Ceverino, Dekel & Bournaud 2010; Krumholz & Burkert 2010). Possibly, a combination of factors is at play. Recent simulations by Krumholz et al. (2017), that include both stellar-feedback and gravitational instability (transport driven) turbulence, show that while both processes contribute to the gas turbulence, transport-driven mechanisms dominate, especially for z > 0.5 and more massive galaxies.

From the observational side, numerous surveys have provided some insights into these hypothesis. Recent work by Johnson et al. (2018) studied the velocity dispersion of ∼450 star-forming galaxies at z ∼ 0.9 from the KROSS survey, finding that the data are equally well explained by a scenario where turbulence is driven by stellar feedback or increased gas fractions. Previous studies have also not provided a definitive answer. In an analysis of z ∼ 3 galaxy analogues (the DYNAMO sample), Green et al. (2014) found a correlation between the SFRs and velocity dispersions that supports the idea that turbulence is driven by stellar feedback. On the other hand, at high redshift (z ∼ 2), the analysis of KMOS^3D survey data by Wisnioski et al. (2015) found only a weak correlation between gas velocity dispersion and SFR or gas fractions at each redshift.

Another possible path to understanding the properties of these young discs is to study their resolved properties at sub-kiloparsec scales. This is particularly challenging when trying to derive intrinsic velocity dispersions. Owing to the limited resolution of observations, the velocity gradient present in rotating discs artificially increases the observed velocity dispersion of galaxies. The effects of beam smearing are particularly problematic in the centres of galaxies, where the velocity field rapidly changes. A few works used adaptive optics to measure the resolved intrinsic velocity dispersion maps of high-z disc galaxies, overcoming some of the issues caused by beam smearing (e.g. Genzel et al. 2011; Swinbank et al. 2012a; Newman et al. 2013). Genzel et al. (2011) find that the velocity dispersion maps are broadly compatible with a constant distribution, as seen in local disc galaxies, despite its higher value (σ ∼ 60  km s⁻¹compared with 10–20 km s⁻¹in local galaxies). However, other studies find a relation between local SFR surface densities (Σ_SFR) and velocity dispersions, which seems to point to some degree of structure (e.g. Swinbank et al. 2012a; Lehnert et al. 2009).

Studying turbulence in the high-z star-forming clumps could also provide some insights into how they differ from their local counterparts. High spatial resolution is necessary to study both the morphology of the velocity dispersion and the clumps’ turbulence. However, many studies of high-redshift galaxies are still severely hampered by the relatively low spatial resolution achieved with the current facilities. A possible strategy to improve this is to target lensed galaxies, where even distant z ∼ 1 galaxies can be resolved down to a few hundred parsecs (e.g. Jones et al. 2010; Livermore et al. 2015; Leethochawalit et al. 2016; Yuan et al. 2017). Here, we target a sample of typical z ∼ 1 clumpy discs that, due to strong gravitational lensing, appear as extremely extended objects in the sky. Their high magnification factors, as well as the high-quality MUSE¹ data acquired with excellent seeing, allows us to alleviate resolution issues, providing a sub-galactic view of the properties of these typical galaxies.

This paper is organized as follows. In Section 2, we present the sample and the MUSE observations used in this work as well as other ancillary data and the lensing models used to recover the intrinsic properties of these galaxies. In Section 3, we derive the integrated physical properties of this sample from the MUSE spectra and HST photometry, and in Section 4 we study the resolved kinematic properties of the galaxies. We conclude with a discussion of the results and comparison with other samples and works in Section 5. Throughout this paper, we adopt a Λ-CDM cosmology with Ω = 0.7, Ω_m = 0.3 and H₀ = 70 km s⁻¹ Mpc⁻¹. Magnitudes are provided in the AB photometric system (Oke 1974). We adopt a solar metallicity of 12 + log (O/H) = 8.69 (Allende Prieto, Lambert & Asplund 2001) and the Chabrier (2003) initial mass function.

2 OBSERVATIONS AND DATA REDUCTION

In this work, we analyse a sample of eight strongly lensed galaxies that lie behind eight different galaxy clusters (see Table 1). All of these clusters were observed with MUSE: Abell 370 (A370), Abell 2390 (A2390), MACSJ0416.12403 (M0416), MACSJ1206.2-0847 (M1206), Abell 2667 (A2667), and Abell 521 (A521) within the MUSE Guaranteed Time Observations (GTO) Lensing Clusters Programme (PI: Richard); Abell S1063/RXJ2248-4431 (AS1063) during MUSE science verification (PI: Caputi & Clément); and MACSJ1149.5+2223 (M1149) during Director's Discretionary Time (PI: Grillo), targeting the newly discovered supernova in a lensed galaxy (Kelly et al. 2015).

We base our selection criteria on the apparent (i.e. due to gravitational lensing magnification) size of the targets with the goal of resolving these galaxies at sub-kiloparsec scales. Within the sample of lensing clusters observed with MUSE, we select the highly magnified (μ > 3) and exceptionally extended (>4 arcsec²) objects, generally dubbed gravitational arcs, that allow to probe properties at physical scales down to a few hundred parsec at z = 1. For comparison, without lensing, at this redshift, a good seeing of 0.6 arcsec results in a resolution of ∼5 kpc. In order to investigate the physical properties of these gravitational arcs, such as kinematics and SFRs, we require that strong, non-resonant emission lines are detected. Therefore, we target objects with visible [O ii ]λ3726, 29 emission in the range 4650 and 9300 Å. This limits our sample to eight galaxies between 0.6 < z < 1.5.

For clarity, we refer to these galaxies by the cluster name plus their respective number in the lensing models used in this work (see Table 1), for example ‘A370-sys1’ or ‘A2667-sys1’. Exceptions to this are the objects in cluster AS1063 and A2390, which do not have multiple images and are referred to as ‘AS1063-arc’ and ‘A2390-arc’.

In the following subsections, we discuss the observations and data reduction and briefly summarize the ancillary data available for these targets, four of which are part of the Frontier Fields HST programme (Lotz et al. 2017). We also present the lensing models and magnification factors used throughout the paper to derive the intrinsic properties of the gravitational arcs.

2.1 Muse observations and data reduction

The targets were observed with MUSE between 2014 and 2017. MUSE is an Integral Field Unit instrument with a field of view of 1 arcmin × 1 arcmin, sampled at 0.2 arcsec pixel⁻¹, and covers the optical wavelength range between 4650 and 9300 Å with a spectral sampling of 1.25Å (Bacon et al. 2014). The GTO targets have a variety of integration depths (between 2 and 6 h). Each observation is comprised of individual 1800 s exposures. The position angle is rotated by 90 deg between each exposure, minimizing the stripe pattern that arises from the image slicer. In order to further minimize this pattern, a small dither (⁠|${<}1{\rm \,arcsec}$|⁠) was added to each pointing. The observing strategy for the non-GTO targets (AS1063, MACS1149, and partially A370) is described in Karman et al. (2015) and Grillo et al. (2016), respectively. They follow very similar strategies to the one used within the GTO sample, including a small dithering and a 90 deg rotation between exposures. A521 was observed in 2017 October with the newly commissioned Adaptive Optics Facility.

The data from all eight clusters was reduced using the ESO MUSE reduction pipeline version 1.2 (Weilbacher, Streicher & Palsa 2016). The calibration files used to perform bias subtraction and flat-fielding, including illumination and twilight exposures, were chosen to be the ones closest to the date of the exposure. To perform flux calibration and telluric correction, all 15 available standard stars were reduced and their flux and telluric response derived using the esorex recipe muse_standard. The response curves were visually inspected and extreme responses (e.g. very high or low, or with unusual features compared with others) were removed before producing a median flux calibration and telluric correction with the six remaining response curves. These median calibration curves were applied to all individual cubes.

The individual cubes were aligned with HST F814W images using a combination of SExtractor (Bertin & Arnouts 1996), to identify bright sources in individual the cubes and HST image, and scamp (Bertin 2006), to calculate the offset between the two.

Sky subtraction was performed using the ESO pipeline and the remaining sky subtraction residuals were removed from each individual exposure using the Zurich Atmosphere Purge tool (zap version 1; Soto et al. 2016), a principal component analysis that isolates and removes sky line residuals. Finally, the individual exposures were combined, rejecting voxels more than 3σ from the median, in order to eliminate cosmic rays. To normalize the exposures, we estimate the sky transparency using the median of the fluxes estimated with SExtractor. For each cluster, the exposure with the highest median flux was taken as the reference and the mean flux of the other exposures was rescaled to match this reference, before combining the data. A second sky residual subtraction using zap was then performed on the combined cubes, and a median filter along the wavelength axis (box of 100 Å) was applied in order to further smooth the background. During both operations, zap and median subtraction, the brightest sources were masked when estimating the background. Finally, Bacon et al. (2014) used the same data reduction and found that the variance propagated by the pipeline is underestimated by a factor of 1.6. We correct the variances in out cubes by this factor.

To determine the point spread function (PSF) size, the final cubes were compared to the available HST data. By weighting the MUSE cubes with the F814W filter transmission curve, we produced MUSE F814W pseudo-broad-band images. We model the PSF following the Bacon et al. (2017) method, which consists of assuming a Moffat profile, with a fixed power index of 2.8, and fitting the full width half-maximum (FWHM) by minimizing the difference between the MUSE pseudo-broad-band and the HST F814W image convolved with the Moffat kernel (F606W was used for A521). With the exception of AS1063, all targets were observed in excellent seeing conditions, always better than a FWHM of 0.70 arcsec (see Table 1).

2.2 Ancillary data and lensing models

The gravitational arcs used in this work were selected to be strongly magnified and are either multiply imaged or very close to the multiple-image region near the cluster core. In order to recover their intrinsic (i.e. corrected for the lensing magnification) properties and morphology, accurate lensing mass models are required.

The selected clusters have all been well studied in the past using numerous and deep HST observations. All clusters except A370, A2390, and A521 are part of the Cluster Lensing And Supernova survey with Hubble (Postman et al. 2012), for which photometry is available in 12 filters (see Table B1 for full list). Furthermore, four of them – AS1063, A370, MACS0416, and MACS1149 – are part of the Frontier Fields initiative (Lotz et al. 2017), that provided deep imaging in the F435W, F606W, F814W, F105W, F125W, F140W, and F160W bands.

Only A2390 and A521 were not part of any of these programmes. For A2390, abundant HST data is also available (see Richard et al. 2008; Olmstead et al. 2014, as well as Table B1). For A521 only WFPC2/F606W band is available plus a NIRC2/Ks band ground-based image retrieved from the Keck archive. This wealth of data has allowed a large number of multiple systems to be identified (between 5 and 60, creating up to 200 multiple images, e.g. Caminha et al. 2017) resulting in very well-constrained lensing models. The MUSE-GTO programme confirmed, as well as discovered, a significant number of these systems, improving the mass models in the process (Richard et al. 2015; Bina et al. 2016; Lagattuta et al. 2017; Mahler et al. 2018).

The mass models used in this work were all constructed using the lenstool² software (Jullo et al. 2007), following the same methodology described in detail in earlier works (e.g. Richard et al. 2009). To summarize, we model the 2D-projected mass distribution of the cluster as a parametric combination of cluster-scale and galaxy-scale dPIE (double pseudo isothermal elliptical; Elíasdóttir et al. 2007) potentials. To limit the number of parameters in the model, the centres and shapes of the galaxy-scale components are tied to the centroid, ellipticity and position angle of cluster members as measured on the HST images. The selection of cluster members is performed using the red sequence in a colour–magnitude diagram (see e.g. Richard et al. 2014) and confirmed with MUSE spectroscopy. Constraints on the lens model are robustly identified multiple images selected based on their morphology, colours and/or spectroscopic redshifts. Although previous lenstool models produced by our team were published for some of the clusters (e.g. Richard et al. 2010b for A2390 and A2667; Richard et al. 2010a for A370), we update them using the most recent spectroscopic information, in particular coming from MUSE (Jauzac et al. 2016; Lagattuta et al. 2017). We summarize the references of the lensing models used in Table 1.

3 SAMPLE CHARACTERIZATION

We start by describing the lensed morphology of these galaxies and then proceed to derive their integrated properties, such as mass, gas-phase metallicity,³ and SFR.

3.1 Morphology

Most galaxies studied in this work have multiple images, i.e. the same object can be seen in two (or more) different positions in the sky with different distortions due to lensing effects (see Fig. 1). Often, the most extended objects do not contain the entire image of the galaxy, that is, only part of the galaxy is lensed into multiple images, which typically appear mirrored with respect to the critical line. One advantage of these multiple images is that the integrated signal of these gravitational arcs is very high, since we observe the same galaxy twice or more, depending on the multiplicity of the image. However, a potential drawback is that properties measured in gravitational arcs, such as mass or morphology, may not be a correct representation of the full galaxy, since the galaxy images may be only partially lensed. To overcome this issue, less magnified but complete images, usually called counter-images, can be used to recover the full morphology of the galaxy and to rescale the physical properties measured in the high signal-to-noise ratio (S/N) data from gravitational arcs.

Open in new tabDownload slide

Open in new tabDownload slide

Figure 1.

Gravitational arcs sample. For each galaxy: No frame: HST composite image with filters F435W, F814W, and F160W in blue, green, and red, respectively (except A521-sys1 where F606W is green and blue, and NIRC2-K band in red). The photometric aperture is plotted in dashed red. Yellow frame: corresponding MUSE [O ii ] pseudo-narrow band of the same sky region. The spectroscopic aperture is plotted in thick blue. Blue frame: F814W band source reconstruction. Both the photometric (HST) and spectroscopic (MUSE) apertures are plotted in red and blue, respectively. Blue frame: F160W band source reconstruction. The magenta cross marks the morphological centre of the galaxy, and the ellipse is placed at 1 effective radius with the inclination and position angle derived from the morphological fit.

A good example of all these effects is A370-sys1, dubbed ‘the dragon’ (see Fig. 1, upper right column). The ‘head’ of the dragon is the complete (or counter) image of this galaxy, where the entire galaxy is lensed with an average magnification factor of 9 and the distortion is small. The ‘body’ of the dragon (the elongated image) is composed of two to four multiple images (Richard et al. 2010a), each of them containing only a part of the galaxy. The magnification in these multiple images is much higher than in the complete image, reaching factors of almost 30. It is in these higher magnification images that the smallest physical scales can be resolved, down to 100 pc for A370-sys1. Similar structures can be seen in A2667-sys1, M1206-sys1, M0416-sys28, and A521-sys1, while AS1063-arc and A2390-arc are single images.

We apply the lensing model to each set of multiple images to obtain reconstructed images, or source plane images images. From these source plane images, we then estimate the effective radius of these galaxies. We do this by fitting a 2D exponential profile to the reconstructed F160W HST image of the complete images. We leave as free parameters the centre of the galaxy, the effective radius, the inclination and position angle. The fit was done using the astropy 2D modelling package (Astropy Collaboration 2013). We caution that these galaxies display complex morphologies, with visible spiral arms and bulges, that were not modelled in the fit. However, we masked the central bright region in A370-sys1, AS1063-arc, M1149-sys1, and A2390-arc, where a bulge is possibly present, to ensure it did not bias the fit. We list the fits parameters in Table 2 and plot ellipses with the best-fitting position angle, inclination and centre at one effective radius with magenta markings in the lower right-hand panels of Fig. 1.

3.2 Spectra extraction and photometry

The rest of the analysis is performed in the image plane, working directly in data cubes and HST images. For the spectral characterization of these galaxies, we are mainly interested in extracting high S/N spectra, where emission lines and continuum can be best analysed. Since the sources have an irregular shape, we do extract spectra from circular or elliptical apertures in the MUSE cube, but instead define extraction areas by imposing a flux threshold on MUSE pseudo-narrow bands. To achieve this, we produce MUSE pseudo-narrow band images centred on the [O ii ] doublet, where we choose a spectral width that maximizes the S/N in the narrow band. We measure the background flux level and variance in these images, and select pixels that are above a 3σ threshold, summing the spectra of the respective spaxels (see Fig. 1). For M0416 and M1149, we exclude multiple images heavily contaminated by cluster members. For M1206, we did not include the counter-image, since it does not significantly increase the S/N.

When coadding the spectra, magnification effects are corrected on a spaxel by spaxel basis, scaling the flux in each spaxel by 1/μ. This is done in order to avoid differential magnification issues. None the less, we note that applying a global (i.e. averaged) magnification factor leads to negligible differences on the overall sample. The average μ value for each spectrum is listed in Table 1. The variance propagated during data reduction is extracted in a similar way.

Photometry was derived using the publicly available HST data. Fluxes are measured over the multiple image which gives a complete coverage of the source in each system. Therefore, for A370-sys1, A2667-sys1, M1206-sys1, and M0416-sys28, the photometry was measured in the counter-image and not in the arc. To minimize contamination from cluster members that are too close to the gravitational arcs, we performed a 2D fit to the light of close cluster galaxies. We do this separately for each photometric band, by modelling each cluster member with a Sérsic profile and subtracting it from the image. The residuals of this subtraction were measured and included in the error budget of the photometry. Next, the F814W image, that matches the MUSE wavelength range, is used to define an extraction mask, by measuring the mean background and variance and placing a 3σ threshold. For A521, the F606W image was used, since there no F814W observations available. We also mask strong residuals arising from the cluster member subtraction. The flux inside the aperture is summed and the background of the image subtracted. The luminosity-weighted magnification factors over the corresponding HST apertures are calculated and the observed fluxes corrected. Uncertainties on the magnification factors are obtained from the Monte Carlo Markov Chain (MCMC) samples produced by lenstool as part of the optimization of the mass model. The photometry and magnification factors of the HST apertures are listed in Table B1 and the spectra are shown in Fig. A1.

The magnification-corrected photometry was used to renormalize the MUSE spectra, correcting both for aperture size and image multiplicity. We do this by measuring the flux in the extracted spectra (applying the F814W transmission curve), and normalizing it to the value obtained from the HST F814W image (discussed below). This process yields high S/N spectra, corrected for lensing and aperture effects, that we use to determine the integrated properties of the sample.

3.3 Emission line measurements

Besides strong [O ii ] λ3726,29 present in all galaxies due to the selection criteria, six out of the eight galaxies display [Ne iii ] λ3869 emission. Other strong lines, such as [O iii ] λ4959,5007 and Balmer emission and absorption lines, as H β and H γ, can also be found in the spectra depending on the object redshift. The narrow line profile of the emission lines, particularly at its centre (see in Section 4), as well as the absence of [Mg ii] emission, make the presence of broad-line AGN unlikely, although we cannot completely rule out this possibility.

In order to better constrain the properties of the emission lines (peak, flux and line width) we first fit and subtract the spectral continuum. To do this, we use the ppxf code version 6.0.2 (Cappellari 2017) and a sample of stellar spectra from the Indo-US library (Valdes et al. 2004), which includes all stars with no wavelength gaps (448 stars). This library covers the entire MUSE observed wavelength range with a spectral resolution of 1.35 Å (FWHM) at rest frame. This corresponds to a FWHM of 2.7 Å in the observed frame at z = 1. This is comparable to the MUSE LSF, about ∼2.5Å, with a wavelength dependence, as determined in Bacon et al. (2017). The continuum fit is performed masking emission lines. To improve the fit, we add a low-order polynomial to the templates and multiply by a first-order polynomial. We manually tune the degree of the additive polynomial by evaluating the quality of the fit as well as the distortion it causes in the shape of the original best fit (i.e. the combination of the stellar templates), since the polynomial can change the original continuum shape. Results of the fit are shown in Fig. A1.

The continuum-subtracted spectra were used to measure emission lines fluxes. We use a Gaussian model implemented in the mpdaf package (Bacon et al. 2016). To estimate the error of the fit, we fit 500 realizations of the spectra, drawn randomly from a Gaussian distribution the with mean and variance corresponding to the observed spectrum flux and variance. We then take the error of the measurements (flux, FWHM and z) as the half-distance between the 16th and 84th percentile of the 500 fits. The flux values are presented in Table B2, and the redshift and line width of the spectra are listed in Table 3.

3.4 Metallicity and star formation rates

From the emission line fluxes, we simultaneously derive the global metallicity and attenuation factor of the galaxies. We do this by fitting multiple metallicity and attenuation sensitive line ratios. This method is a modification to that presented by Maiolino et al. (2008, hereafter M08), whereby we extend their line ratio list and place it on a more formal Bayesian footing. We use the following line ratios:

1. |$[\rm{O\,{\small II}}\,]\lambda \,3727$|/H β, M08

2. |$[\rm{O\,{\small III}}\,]\lambda \,5007$|/H β, M08

3. |$[\rm{O\,{\small III}}\,]\lambda \,5007$|/[O ii ]_λ 3727, M08

4. [Ne iii ]/[O ii ] λ 3727, M08

5. R23⁴, M08

6. [O iii ] λ 5007/[O iii ] λ 4959, 2.98 (Storey & Zeippen 2000)

7. H γ/H β, 0.466 (Osterbrock & Ferland 2006, case B)

8. H δ/H β, 0.256 (Osterbrock & Ferland 2006, case B)

9. H δ/H γ, 0.549 (Osterbrock & Ferland 2006, case B)

10. H7/H γ, 0.339 (Osterbrock & Ferland 2006, case B)

The line ratios [1–5] are the same as those adopted by M08. These line ratios provide constraints on metallicity. But, because line ratios [1,3,5] compare lines at disparate wavelengths, they are also sensitive to the attenuation. So, to help alleviate the metallicity-dust degeneracy, we include ratios that are predominantly independent of metallicity [6–10].

We use all line ratios for which the lines are observed. However, since H β is not available for galaxies at z >0.9, we replace the [O ii ]/H β in ratio [1] with [O ii ]/H γ (assuming case B H γ/H β ratio of 0.466, for T_e = 10 000 K and low electron density). In Table 4, we list the specific line ratios used for each galaxy.

Figure 2.

Example of simultaneous determination of metallicity and extinction for galaxy AS1063-arc. The marginalized distributions for metallicity and extinction are displayed in the top left and bottom right panel, respectively. The bottom left panel displays the joint distribution. Similar plots for the remaining galaxies can be found in Appendix D.

Open in new tabDownload slide

M1149-sys1 was left out of this analysis, since only [O ii ] is present in the MUSE spectrum. All the galaxies analysed have well determined metallicities and attenuations, ranging from 8.71 to 9.05 and 0.09 to 2.86, respectively, with symmetrical marginalized distributions (except for attenuation in M0416-sys28, A2390-arc, and A521-sys1).

For A2390, the only metallicity diagnostic that could be used was |$[\rm{O\,{\small II}}\,]\lambda \,3727$|/H γ, since the S/N does not allow a robust measurement of the [Ne iii ] flux. As the |$[\rm{O\,{\small II}}\,]\lambda \,3727$|/H γ diagnostic is degenerate, with both low- and high-metallicity solution, we investigate each branch by repeating the fit allowing only low or only high metallicity values ([7.0,8.6] and [8.6,10.0] priors, respectively). For solutions found in each branch, the [Ne iii ] flux would be ∼6.2 × 10⁻¹⁸ and ∼1.2 × 10⁻¹⁸ erg cm⁻² s⁻¹ for the low- and high-metallicity solutions, respectively. The low-metallicity estimate approximately corresponds, for example, to the flux of H δ seen in this galaxy (∼5.7 × 10⁻¹⁸ erg cm⁻² s⁻¹), so the low-metallicity branch seems less likely than the high-metallicity solution. The high-metallicity solution would yield lower [Ne iii] fluxes that would likely be undetected given the ∼0.3 × 10⁻¹⁸ erg cm⁻² s⁻¹ 1σ spectral noise for this galaxy. Furthermore, the SFRs derived via the Balmer lines and the [O ii ] line (corrected for metallicity) are closer to the high-metallicity case than the low one. These indications do not rule out the low-metallicity hypothesis, but we adopt the high-metallicity values as the most probable solution for A2390-arc throughout the rest of this work.

The SFR are calculated both from Balmer lines and the [O ii ] doublet. For the Balmer lines, the Kennicutt (1998a) calibration is adopted, adapting this calibration to be used with H β or H γ and assuming their ratios relative to H α are given by Case B theory at electronic temperature T_e = 10 000 K and low electron density. For the [O ii ], the metallicity-dependent Kewley, Geller & Jansen (2004) calibration is used. The SFR from Balmer lines is obtained by drawing 500 random values of H β or Hγ fluxes (assuming a Gaussian distribution) and correcting them for dust attenuation using a sample of 500 τ_v values drawn from the marginalized distribution obtained from the metallicity and attenuation fit. The dust corrected line fluxes are converted to SFR using the Kennicutt (1998a) calibration and the SFR taken as the median of the sample. For the SFR estimated from [O ii ], an equivalent process is used. We produce a sample of [O ii ] fluxes, τ_v and metallicity values and apply the Kewley et al. (2004) calibrations. Finally, both SFR distributions are convolved with the magnification error, since they depend on the absolute flux of the emission lines. We list the calculated SFR in Table 4.

3.5 Stellar mass

We fit both spectra and photometry using the FSPS (flexible stellar population synthesis; Conroy, Gunn & White 2009) models implemented within a Bayesian framework in the prospector⁵ code (Johnson & Leja 2017). The only exception is A521-sys1, for which only the spectrum was fit and used to normalize the spectra, since only one HST filter was available. prospector allows us to fit models produced with FSPS ‘on the fly’, i.e. the fit are not based on a precomputed grid of models but are instead computed for each exact set of parameters tested. We assume a decaying star formation history and choose the Padova isochrones (Marigo & Girardi 2007; Marigo et al. 2008) with the Chabrier (2003) initial mass function and the MILES stellar library (Sánchez-Blázquez et al. 2007; Falcón-Barroso et al. 2011) to generate the models for all fits. We also adopt the Charlot & Fall (2000) dust attenuation here, with a fixed index of −1.3 for the star-forming regions, applied only to the star-light originating in stars younger than 10 Myr, and an index of −0.7 for the global dust screen as suggested by Charlot & Fall (2000), applied to all starlight equally.

In total, six physical parameters were fit: stellar mass, the e-folding time of the star formation history, stellar population age, stellar metallicity and optical depth of both birth clouds and ISM (see Table 5 for details and priors). Beside these physical parameters, a third-order polynomial is also fitted (similarly to what is done with ppxf). For the sake of simplicity on an already quite complex fit, we decided to mask the emission lines during these fits, and do not attempt to fit them during the continuum fit with prospector and do not include nebular emission in the models.

The fits for each galaxy can be found in Appendix C. The best values and errors, estimated from the marginalized distributions as in the previous section, can be found in Table 5. The magnification errors were included at this point in the mass and SFR error estimates. Overall, the best fits provide a good match to both the spectral and photometric data. Most of the derived model parameters are well constrained, with the exception of the attenuation, which in some instances have large tails extending to high values. In principle, the value of τ_BC obtained here should agree with the attenuation derived from the emission lines. In practice, since no FUV data is available, nor were the emission lines fitted with prospector, this value cannot be constrained and will be degenerate with τ_ISM. The obtained stellar masses vary from ∼1 ×  10¹⁰ to ∼9 ×  10¹⁰ M_⊙, with a mean value of 3.65 ×  10¹⁰ M_⊙, and with an e-folding time that varies from 3 to 9 Gyr and dominant ages from 3 to 7 Gyr.

3.6 Comparison with other samples

We can now use the physical parameters previously derived both from emission lines and photometry to place this sample in context. For M1149-sys1, only the [O ii ] line is present in the MUSE wavelength range, so it is not possible to constrain metallicity and SFR from these data. Two metallicity estimates exist for this galaxy in the literature. Yuan et al. (2011) derive a metallicity of 12 + log (O/H) = 8.36 ± 0.04 using the [N ii ]/H α ratio, while Wang et al. (2017) measure 12 + log (O/H) = 8.70 ± 0.10 using HST grism spectra. Wang et al. (2017) use M08 calibrations, including line ratios with [O ii ], [O iii ], and H β emission lines, and calculate metallicity using a Bayesian inference framework, comparable to what is done here. We make use of both their metallicity and SFR (16.99 ± 4.3 M_⊙ yr⁻¹) values derived from emission lines throughout this work.

First, we compare the sample with the mass–SFR relation (top panel of Fig. 3). Within errors, all the galaxies lie close to what is measured in main-sequence star-forming galaxies. The only exception is A370-sys1 that has a slightly lower SFR than expected from its derived mass. We then compare this sample with the fundamental metallicity relation (FMR; Lara-López et al. 2010; Mannucci et al. 2010), that defines a tight relation between stellar mass, gas metallicity, and SFR of normal star-forming galaxies (bottom panel of Fig. 3). Five out of the seven arcs – AS1063-arc, A370-sys1, A2390-arc, M0416-sys28, and M1206-sys1 – differ by up to 0.1 dex from this fundamental relation, which is higher than the dispersion derived for the local SDSS sample of about 0.05 dex, but still within the 0.5 dex dispersion observed for higher redshift samples (Mannucci et al. 2010). We conclude that all the galaxies studied in this work are normal star-forming galaxies and are good representatives of the population of disc galaxies at z ∼ 1.

Figure 3.

Top panel: Gravitational arcs sample plotted against the main sequence. We plot in solid grey lines the parametrization of Whitaker et al. (2012) for this relation at z = 0.6 and z = 1, with the dispersion in grey. Bottom panel: Distance to the FMR from Mannucci et al. 2010 (with the 1σ range of the SDSS sample form the same work shown in grey). Most galaxies lie close to this fundamental relation. The same colour scheme is used in both plots to identify objects. Both relations considered, these galaxies are typical star-forming galaxies for z ∼1.

Open in new tabDownload slide

We can estimate the gas mass of these galaxies by using the Kennicutt–Schmidt (KS) law (Kennicutt 1998b) to convert SFR surface density to gas surface density. We take the global (dust and magnification corrected) SFR derived from the Balmer lines of the integrated spectrum and we calculate the SFR density by dividing it by the photometry aperture area (the area to which the spectra were normalized). Note that the area is calculated in the source plane. The SFR is then converted to gas surface density, and thereby also the total gas mass. These are listed in Table 6. The SFR densities, gas masses (M_{g,  KS}), and gas fractions (M_{g,  KS}/(M_{g,  KS} + M_⋆)) are listed in Table 6. These galaxies have low gas fractions that range from ∼0.1 to 0.3, despite the star-forming clumps clearly visible in most of them (particularly AS1063-arc, A370-sys1, and M1206-sys1). These values are consistent with the values obtained from the DYNAMO survey, a sample of comparable local galaxies that covers a stellar masses of 10⁹–10¹¹ M_⊙ and SFRs of 0.2–100 M_⊙ yr⁻¹, selected to be analogues to z ∼ 3 galaxies (Green et al. 2014). The gas fraction in our galaxies are also similar to what White et al. (2017) found from making CO measurements and inverting the KS law.

4 RESOLVED PROPERTIES

4.1 Observed velocity maps

The velocity maps (or kinematic maps) were derived measuring the [O ii ] emission lines present in all galaxies. The fits were performed with a slightly modified version of the camel⁶ code (Epinat et al. 2010), where we added the option of fitting binned data. camel fits emission lines with 1D Gaussian models, having as free parameters the redshift, the Gaussian FWHM and the flux of each emission line. In order to more robustly probe the outer parts of the galaxies, we bin the data using the Voronoi binning method⁷ presented in Cappellari & Copin (2003). We bin on the S/N of the [O ii ] pseudo-narrow bands used in Section 3.2. We impose a S/N of 5 per bin in all galaxies, which we find to be a good compromise between high spatial sampling and robust spectral fits. The velocity dispersion was allowed to vary between 0 and 250 km s⁻¹. The 2D observed velocity maps of the eight galaxies are presented in Fig. 4. These maps were inspected and the velocities values were corrected (adding or subtracting a constant value in the entire field) so that the velocity at the morphological centre (as given by the 2D exponential fit to the source plane image, see Table 2) corresponds to the kinematic centre and has velocity zero.

Open in new tabDownload slide

Open in new tabDownload slide

Figure 4.

Observed and model velocity maps of the eight gravitational arcs. For each galaxy, the three first upper panels display the velocity maps: the observed velocity derived from [O ii ] emission, the best model (with the model name) and the residuals. On the first panel, a black circle gives the seeing size for that observation. The three lower panels are the equivalent for the velocity dispersion. Both colour bars are in units of  km s⁻¹. We have rejected and masked (in grey) values with unreliable beam-smearing corrections (σ_smear > σ_obs)

4.2 Fitting velocity maps

We fit the 2D velocity maps with the three following kinematic models: an arctangent model (Courteau 1997), an isothermal sphere model (Spano et al. 2008) and an exponential disc model (Freeman 1970). The first is an empirical model commonly used to fit local galaxies. The isothermal sphere model is a good approximation of the kinematics of a dark matter dominated galaxy, under the assumption that the dark matter halo has an isothermal profile. The exponential disc assumes a smooth exponentially decaying distribution of baryonic matter (and no dark matter) as seen in most disc galaxies (see Epinat et al. 2010 for a summary of these models and analytical expression implemented in this work).

We fit these models directly to the 2D velocity maps produced by camel (i.e. in image plane). Note that we do not fit the velocity dispersion maps. To account for lensing distortions, we produce high-resolution displacement maps using the optimized cluster mass models and the lenstool software. These maps describe the geometrical transformation from source to image plane by predicting the spatial shifts (or displacement) that have to be applied to each source plane pixel to be projected to the image plane image.

The first step of the fit is to create a 2D kinematic map in source (unlensed) plane. We apply the lensing distortion to this kinematic model, using the displacement maps. We then convolve the lensed model velocity map with the seeing, using a Moffat kernel, flux-weighted by the [O ii ] flux image from camel. Finally, we apply the same binning to the models as the observations before comparing the modelled and the observed velocity fields. It is this final product that is compared with the data.

Each model has seven free parameters, that are explored using a MCMC sampler:

• X,Y position of the kinematic centre in image plane.

• V_sys systemic velocity.

• Inclination of the galaxy (inc).

• Major axis kinematic angle in source plane (PA).

• Maximum velocity (V_max).

• Transition radius (r_t).

The first three parameters (X,Y,V_sys) are mainly introduced to correct for small inaccuracies in the lensing correction (V_sys is used to keep the v = 0  km s⁻¹ at the centre even if the centre changes slightly). We assume that the kinematic centre is the same as the morphology centre (measured in the reconstructed image), and correct the observed velocity fields in order to have zero velocity at this position. To correct any possible errors when passing from the position in source plane to image plane, during the fit we allow the kinematic centre to vary within a 2×2 MUSE pixel box around the original position and let the velocity adjust up to ±100  km s⁻¹. The corrections for all galaxies were small with displacements of 1 pixel and velocity adjustments of up to ±30  km s⁻¹.

The maximum velocity and transition radius (V_max, r_t) correspond to the maximum velocity (and the radius at which it occurs) for the exponential disc and isothermal sphere models. The arctangent model velocity does not reach a maximum at a finite radius, so the transition radius in the arctangent model corresponds to the radius at which 70 per cent of the asymptotic velocity (V_max) value is reached, and it is therefore generally lower for the arctangent model than for the other two kinematic models.

We choose a wide prior for the maximum velocity and transition radius that is the same for all galaxies. The inclination was given a wide prior of ±20 deg around the inclination derived from the morphological fit. The only exception is AS1063-arc, observed with the worst seeing in this sample (∼1 arcsec), where the inclination was unconstrained by the kinematic data and was set to the value measured in the HST image. The position angle was also given a wide prior (±20 deg) around the position angle measured in the exponential fit to the HST images. Finally, we present the Bayesian evidence of each model, a relative measure of how well a model performs, as calculated by emcee.

4.3 Results

The results for each model are presented in Table 7, where the best value of each parameter corresponds to the median of the marginalized distributions. Also here, the half distance between the 16th and 84th percentiles was taken as the error for each parameter. From the Bayesian evidences values, the preferred kinematic model in this sample is the arctangent model, with galaxies AS1063-arc, A2667-sys1, M0416-sys28, and M1206-sys1 being better fit by this model. (However, notice that the fit of M0416-sys28 is not constrained, with the transition radius hitting the upper limit allowed.) The remaining galaxies, A370-sys1, A521-sys1, and M1149-sys1 are better fit by an exponential disc model and only A2390-arc by an isothermal sphere. And although there is a clearly preferred model for each galaxy, none of the models is globally preferred. The position angles (PA), that are expected to be the same among the three models, agree quite well for each galaxy. The inclination shows a larger spread, but still displays a reasonable agreement. The 2D maps for the best model with the parameters of Table 7 are plotted in the middle panel of Fig. 4.

As previously stated, intrinsic velocity dispersions are difficult to estimate due to the beam-smearing effect. Here, we benefit from high magnification of the lensing and the good seeing of the MUSE observations that reduce this problem, but nevertheless we still correct each pixel using the 2D velocity kinematics. From the best model for each galaxy, we calculate the expected velocity dispersion created by the beam-smearing effect (σ_smear) applying equation A24 from Epinat et al. (2012). At the [O ii ] wavelengths for this sample, MUSE has a line spread function (LSF) with FWHM of ∼40–80 km s⁻¹ (as from Bacon et al. 2017 fit of the UDF mosaic data). We subtract in quadrature the LSF value at the [O ii ] wavelength of each galaxy. The intrinsic velocity dispersion, corrected for instrument broadening and velocity smear, is then given by |$\sigma = \sqrt{\sigma _{\textrm{obs}}^2 - \sigma _{\textrm{inst}}^2 - \sigma _{\textrm{smear}}^2 }$|⁠. We plot these 2D ’ panel in Fig. 4.

All galaxies are rotation-dominated, with regular 2D velocities and V/σ ratios between ∼2 and ∼10, estimated from the maximum velocity V_max and the mean of the σ maps (see Table 8). These values are in agreement with V/σ ratios found for other high-z samples such as the KMOS^3D kinematics survey (Wisnioski et al. 2015), with V/σ between ∼3 and 11 for galaxies at z ∼ 1, the MUSE sample of low-mass galaxies (Contini et al. 2016), with the range ∼1–5 and the KLASS survey of lensed galaxies at 0.7 < z < 2.3 (Mason et al. 2017), with V/σ ∼8 at z ∼ 1.

5 DISCUSSION

In this section, we discuss the implications of the kinematic fit results and the morphology of the velocity dispersion maps. We also explore the connection between star formation and turbulence in these typical z ∼ 1 discs. A summary of the individual properties of the galaxies, and comparison with previous studies of these objects, can be found in Appendix E.

5.1 Best kinematic model

Recent work by Genzel et al. (2017) and Lang et al. (2017) found that both the individual velocity curves of massive galaxies and the stacked velocity curve of 101 star-forming galaxies in the range 0.6 < z < 2.6, display a significant drop in the velocity curve beyond the turnover radius. This would imply a higher baryonic mass fraction in high-z galaxies than previously though. We briefly investigate this issue.

As stated before, four out of eight galaxies are better fit by the arctangent model, an empirical model. The exponential disc model, which describes the rotation velocity of a gas disc, has lower Bayesian evidences on three of the galaxies and the isothermal sphere model, an approximation to what is expected for an isothermal dark matter halo, on only one. Most importantly, all models are remarkably similar in the inner parts of the galaxies. For most galaxies, the three models (with the best parameter estimations) are visually indistinguishable up to 1 R_e (e.g. A370-sys1 and A2667-sys1), or even 2 R_e (e.g. AS1063-arc).

In none of the galaxies analysed here do we see a drop in the velocity curve at large radii, which is confirmed by the fact that the arctangent model (the only one that does not decrease velocity with radius) is the most preferred model. However, with the exception of AS1063-arc and A2667-sys1, our data does not extend out as far as the stacked sample of Lang et al. (2017), so we cannot rule out the possibility of a drop in the velocity beyond two effective radii. We also caution the reader that our analysis is too simplistic to explore this subject, since a combination of the two components – baryons and dark matter – is expected, and a joint fit should be done in order to gain some more insight on this matter, as done in Lang et al. (2017).

5.2 Velocity dispersion morphology

We now turn our attention to the velocity dispersion. The mean intrinsic velocity dispersion of these disc galaxies ranges from 15 to 84  km s⁻¹, with a mean value of 48 km s⁻¹. We plot these values in Fig. 5 and compare them with other samples. The mean velocity dispersions are comparable to the global values derived for the KMOS^3D, the MUSE samples for 0.2 < z < 1.4 galaxies, and the lensed sample from the KLASS survey (Wisnioski et al. 2015, Contini et al. 2016, and Mason et al. 2017, respectively). Our velocity dispersions are also in agreement with the values derived for the DYNAMO sample of local analogues (Green et al. 2014; Oliva-Altamirano et al. 2018). When compared to the GHASP survey (Epinat et al. 2010), a sample of nearby late-type galaxies, these lensed galaxies have velocity dispersions higher by a factor of 1.4–4.2.

Figure 5.

Comparison of the velocity dispersion of the gravitational arcs with other samples. Velocity dispersion against stellar mass on the top panel and velocity dispersion against star formation surface density on the bottom (M1149-sys1 is omitted from this plot since we do not derive its SFR). The mean value for the local sample of rotating galaxies from Epinat, Amram & Marcelin (2008) in plotted as a black line, with scatter masked by the grey area; the local (z = 0.1) DYNAMO sample of Green et al. (2014) in black squares (we plot only the rotating disc galaxies); the z ∼ 0.2–1.4 sample from Contini et al. (2016) in black triangles; a sub-set of the z = 1 KMOS^3D galaxies, with velocity dispersions derived by Di Teodoro, Fraternali & Miller (2016) as black crosses (only in the top panel); the lensed galaxies from the KLASS sample (Mason et al. 2017) are plotted as black diamonds (we have excluded possible mergers and galaxies with only upper limits). The KMOS Redshift One Spectroscopic Survey Sample at z ∼ 1 from Johnson et al. (2018) is also plotted in grey dots. The gravitational arcs sample velocity dispersions are in good agreement with all these samples.

Open in new tabDownload slide

Several studies have investigated the kinematics of lensed galaxies (e.g Jones et al. 2010; Livermore et al. 2015; Leethochawalit et al. 2016; Mason et al. 2017), mostly focusing on the velocity fields of these galaxies and not so much in the resolved velocity dispersions. The sample presented in the current work was selected to be extremely extended and with high physical resolution and the high S/N of our data results in relatively small formal errors per bin in the 2D velocity dispersion maps. These are derived with camel and, for AS1063-arc, A370-sys1, A2390-arc, M0416-sys28, A2667-sys1, and M1206-sys1, have a median value that varies between 5 and 8  km s⁻¹(but can be as high as 50  km s⁻¹). The remaining galaxies, A521-sys1 and M1149-sys1, have higher median errors, 24 and 16  km s⁻¹, respectively. Another source of uncertainty is the LSF correction, that we perform using expression 7 of Bacon et al. (2017). For our lowest redshift galaxy, the LSF dispersion is about 1.1Å (see fig. 15 of Bacon et al. 2017), which translates in a dispersion of ∼50 km s⁻¹. For the higher redshift galaxies this value is lower. Neither of these sources of uncertainty directly depends on the position of the bin in the velocity map. On the other hand, the beam-smearing correction does depend on the position of the bin in the 2D map. Calculate the errors associated with this correction is computationally expensive (since a reasonable number of velocity dispersion maps would have to be produced), so we estimate the uncertainties associated with the beam-smearing correction in a simplified way. We produce 16 velocity dispersion models with combinations of perturbed V_max, r_t, PA, and the inclination. Since this is a simplified estimation of the error, we perturbed these parameters by ±20  km s⁻¹, ±1 kpc, and ± 5 deg, values that are considerably higher than our typical errors. Comparing these perturbed models with the best model for each galaxy, we obtain median differences between 5 and 10  km s⁻¹. Overall, our major source of error is the LSF estimation. Nevertheless, Bacon et al. (2017) measure the LSF to be largely constant in the MUSE FOV, with an average standard deviation of ∼0.05 Å, therefore the effects on the study of the velocity dispersion structure are expected to be small.

The observed, σ_smear and intrinsic velocity dispersion 2D maps of the eight galaxies analysed in this work can be seen in the lower panels of Fig. 4. Both the observed and the σ_smear (the dispersion expected due to observational effects) maps have a higher velocity dispersion at the centre of the galaxies (exceptions are A2667-sys1 and M0416-sys28). The intrinsic 2D velocity dispersion maps reveal no significant 2D structure at our precision limit, with the dispersion of the σ maps varying between 7 and 33  km s⁻¹(see Fig. 5). Our intrinsic velocity dispersion maps are broadly consistent with a constant velocity dispersion value throughout the entire galaxy.

In order to better quantify how homogeneous the velocity dispersion maps are we also measure their Gini coefficient (Gini 1912). This is a measure of inequality in a data set and can vary from 1, absolute inequality, to 0, absolute equality (where, in this case, all the pixels have the same value). We obtain values between 0.14 and 0.30 with a mean value of 0.23 (see Table 8). This suggests a fairly homogeneous distribution of the values of the velocity dispersion and it is significantly lower than the values obtained doing the same exercise with the [O ii ] pseudo-narrow bands, that map the star-forming clumps, where we measure Gini coefficients between 0.82 and 0.90. The relatively homogeneous velocity dispersion maps, compared with clumpy [O ii ] and UV morphologies, suggest that the star-forming clumps do not strongly disturb the local velocity dispersion, i.e. the clumps do not seem to locally increase turbulence. We caution however that we did this exercise based on individual pixel values, not deriving velocity dispersions per clump. Recent work analysing the DYNAMO sample by Oliva-Altamirano et al. (2018) found statistically significant structure in six of these local analogue galaxies with increased velocity dispersion possibly associated with clumps. These structures generally show an increase of the order of ∼20 km s⁻¹. This is very close to the precision we can reach with our data, and it is unlikely that we would detect this, even if this structure were present.

5.3 Turbulence and star formation

As stated in the introduction, the reasons why high-redshift discs are more clumpy and turbulent, and how these two properties relate to each other, are still under discussion. While some works point to the higher gas fractions of these younger discs as the reason for their clumpy morphologies and turbulence (e.g. Genzel et al. 2011; Livermore et al. 2015), others found a correlation between SFR and higher turbulence (e.g. Swinbank et al. 2009; Green et al. 2010). We investigate this in the two lowest redshift galaxies of this sample, where the resolved extinction and extinction corrected SFR maps can be robustly derived. We will present and discuss the metallicity maps of these galaxies in a separate paper (Patrício et al. in preparation) and focus here on the SFR density.

To produce the star formation density (Σ_SFR) dust corrected maps, we first subtract the continuum from the observed cube, since this could bias the measurement of the Balmer lines. This is done by fitting the local continuum with ppxf. The data were spatially binned using H γ pseudo-narrow bands, which guarantees that both H γ and H β can be reliably measured in each bin. The line fluxes are fit in the continuum subtracted cube using camel, as done for the observed velocity maps. We calculate the extinction and SFR from the emission lines in each bin as done in the integrated spectrum (see Section 3.4), with the exception that only the |$[\rm{O\,{\small III}}\,]\lambda \,5007$|/H β, |$[\rm{O\,{\small III}}\,]\lambda \,5007$|/|$[\rm{O\,{\small II}}\,]\lambda \,3727$|⁠, |$[\rm{O\,{\small III}}\,]\lambda \,5007$|/|$[\rm{O\,{\small III}}\,]\lambda \,4959,$| and H γ/H β ratios were used and that the parameter space explored this time was smaller (extinction between 0 and 2.5) based on the results obtained in the integrated spectra. The resulting maps can be seen in Fig. 7. Since velocity dispersion and SFR density were derived with a different binning, in the following analysis we only consider the binning in Σ_SFR and obtain the corresponding velocity dispersion by binning the intrinsic velocity dispersion maps to the same grid.

There are different results in the literature on the relation between SFR density and velocity dispersion. While Genzel et al. (2011) found little correlation between the two in their sample of massive giant clumps at z ∼ 2 (Swinbank et al. 2012b), using a sample of galaxies observed with Adaptive Optics, found it was well described by the relation |$\sigma \propto \Sigma _{\rm SFR}^{1/n} + k$|⁠. Since clump selection depends heavily on the criteria and data used (e.g. Dessauges-Zavadsky et al. 2017) we use the full 2D information of these two galaxies and plot in Fig. 8, for each Voronoi bin, the star formation density and the respective velocity dispersion. We also gather data on the velocity dispersion and star formation density of the high-z clumps available in the literature: five giant star-forming clumps at z ∼ 2 from Genzel et al. (2011), from the z ∼ 1 lensed galaxies clumps from Livermore et al. (2015) and from the clumps of nine 0.8 < z < 2.23 galaxies from the HiZELS survey (Swinbank et al. 2012b). for the last two, the Σ_SFR were computed from the SFR and the radii of the clumps.

We first investigate if these two quantities, Σ_SFR and σ are correlated by computing the Spearman rank correlation (ρ) between the two (see Table 9). Since the data on the star-forming clumps gathered in the literature can be quite different from the individual bins studied here, we do this analysis separately.

The full data set on properties of clumps from Genzel et al. (2011), Swinbank et al. (2012b), and Livermore et al. (2015) shows very little correlation between Σ_SFR and σ, with ρ = −0.12. Taking only the data of AS1063-arc, we find a Spearman correlation of −0.21 and for A370-sys1 a correlation of 0.20, once more both very weak.

Also here we analyse the data separately, fitting data only from individual star-forming clumps from the literature, only from AS1063-arc or A370-sys1, and finally all of it. The fit results and their reduced χ² are listed in Table 9. We generally find very low γ values, all below 0.1 except A370-sys1 (0.307 ± 0.041), far from the ∼0.74 found by Swinbank et al. (2012b) and implied by the local KS law. Lehnert et al. (2009) also present pixel-per-pixel data on 11 z ∼ 2 star-forming galaxies observed with AO and find that γ = 0.33 is a good match to their data, which is in agreement with what we see for A370-sys1, but not to AS1063-arc or the individual clumps.

Because we made use of fully resolved (as opposed to clump) data of highly magnified objects, we reach lower SFR surface densities. Indeed, there is only a small overlap in Σ_SFR between AS1063-arc and A370-sys1 when compared with the clump sample gathered from the literature. Other resolved studies, such as Lehnert et al. (2009) and Swinbank et al. (2009), generally reach SFR surface densities as low as 0.1 M_⊙ yr⁻¹ kpc⁻², while much of the our data sits below this limit (although see Cresci et al. 2009 for other Σ_SFR < 0.1 M_⊙ yr⁻¹ kpc⁻² measurements of high-z galaxies). However, for Σ_SFR < 0.1 M_⊙ yr⁻¹ kpc⁻², velocity dispersions are expected to almost not vary with star formation surface density, even in a stellar-feedback driven turbulence scenario (Krumholz et al. 2017), so it is not surprising not to find a correlation at low SFRs. We can nevertheless look at the expected velocity dispersions for different scenarios for the measured stellar formation rate densities. For Σ_SFR between 0.02 and 1 M_⊙ yr⁻¹ kpc⁻², where our data is more reliable, we find values of σ of 20–100  km s⁻¹. Comparing with theoretical predictions from Krumholz et al. (2017, their fig. 4), we see that these high-velocity dispersions for low star formation density are only reached for models that include transport (i.e. gravitational instability) driven turbulence. Our data best compares with the transport+feedback and no-feedback predictions for local dwarf galaxies. These predictions are however done for integrated properties of the galaxies, and not resolved as presented in this work, which makes a direct comparison difficult.

We try to extend this study beyond these two galaxies, where we can correct each bin for extinction, to the other galaxies of the sample, but in a more simplistic way. Since we cannot derive resolved extinction maps for the remaining six objects, we test only if there is a correlation between the observed [O ii ] surface brightness (that probes star formation and is unaffected by magnification since gravitational lensing conserves surface brightness) and the velocity dispersion. A similar analysis, using H α surface brightnesses, was presented in Lehnert et al. (2009) that found a positive correlation between the two. We plot these two quantities and the Spearman coefficients in Fig. 9. To make sure that the Spearman coefficients were not biased by unreliable measurements we have removed data points for which (i) formal errors were 0  km s⁻¹; (ii) in which the velocity dispersion reached the maximum allowed value of 250  km s⁻¹; or (iii) that have S/N lower than 2. These unreliable values are plotted in black squares in Fig. 9.

This crude approach has limitations. Comparing the ρ coefficients obtained for AS1063-arc and A370 using both methods, we see they are extremely different. However, in neither analysis do we find a strong correlation between SFRs (or [O ii ] surface brightness) and turbulence. All the seven considered galaxies have ρ ≤ 0.26, denoting a very faint dependence between these two quantities (M0416-sys28 does not have enough data to perform this exercise). Bearing in mind the limitations of this test, we conclude that there is no strong evidence of local velocity dispersions being sustained by (low) star formation, which is corroborated by the more careful analysis of AS1063-arc and A370-sys1, and is in agreement with previous results from Genzel et al. (2011) and Livermore et al. (2015), but in contrast with what was found by Lehnert et al. (2009), Swinbank et al. (2012b), and Lehnert et al. (2013).

5.4 Energy budget

We calculate both these energies within the photometric aperture (to which the total spectra were normalized), making use of the gas densities from Table 6. The velocity in expression (6) was calculated for the photometric aperture radius using the best model and parameters for each galaxy. We warn that Σ_g and SFR are not independent, since Σ_g was derived using the inverse KS law. We report these values in Table 10.

We find that in our sample, star formation can provide the necessary energy to sustain turbulence, requiring only 6–46 per cent of the available energy. This is somehow at odds with what is found by Krumholz et al. (2017) in their disc evolution model where at z > 0.5 disc turbulence is dominated by transport-driven mechanisms (e.g. from gas accretion) and not stellar-feedback, especially for more massive galaxies, although both processes co-exist. Such a difference can come from the approximations made to calculate these energies, particularly in estimating the gas surface density.

6 SUMMARY AND CONCLUSIONS

We discussed the physical properties of seven clumpy discs at 0.6 < z < 1.5, that are highly magnified by strong lensing and can therefore be resolved to sub-kiloparsec scales. From high-quality MUSE data (with <0.7 arcsec seeing) and abundant HST archival data we derive the following:

• We derive stellar masses, global metallicities and extinctions from multiband HST photometry and multiple line diagnostics for the full sample and for the first time for three of the eight objects (AS1063-arc, A370-sys1, and A521-sys1). Comparing these objects with the main-sequence of star-forming galaxies, we confirm that these are typical star-forming galaxies for the considered redshifts (Fig. 3).

• The observed velocity fields, derived from the MUSE data cubes, all reveal smooth velocity gradients with a velocity dispersion peak at the centre of rotation, typical of rotating discs. These are rotation-dominated galaxies with V/σ between ∼2 and ∼10, also typical of what was found by IFU surveys (Fig. 5). We do not see any evidence of a drop in the velocity curve within the sample, even for galaxies where two effective radii are reached (Fig. 6).

• The observed velocity maps were fit with three different kinematic models (arctangent, isothermal sphere and exponential disc model), with a novel method that allows to fit several multiple images at the same time, accounting for lensing distortions. The arctangent model was a better fit to four galaxies of the sample, the exponential disc model to three, and the isothermal sphere model to the remaining galaxy. However, all models deliver very similar results, particularly in the central part of the galaxies, which seems to indicate that the velocity is not probed far enough to reliably infer which model is the best (and therefore to verify if the velocity curves decay at larger radii).

• The intrinsic velocity dispersion maps (that we correct for beam-smearing effect) do not display any strong 2D structure, revealing fairly isotropic 2D maps with a mean Gini coefficient of 0.23 (compared with 0.87 for the [O ii ] maps, Table 8 ).

• There is no strong evidence of increased (or decreased) velocity dispersion at the positions of the star-forming clumps. We also do not find a strong correlation between the local (resolved) SFR density and the velocity dispersion for the two galaxies where the dust attenuation corrected SFR maps could be derived in a resolved way (Fig. 8). Comparing the resolved results with theoretical predictions, we conclude that the high-velocity dispersion found for relatively low SFR surface densities are more in agreement with scenarios where the turbulence is driven by gravitational instability (Krumholz et al. 2017).

• However, from an energetic point of view, the observed star-formation is sufficient to explain the observed turbulence.

• We find that the star-forming clumps share the dynamical properties of the ionized gas surrounding them (rotation velocity and turbulence). In a future work, we will study the chemical composition of these galaxies and the clumps within them, to further test the in situ clump formation hypothesis.

Figure 6.

Radial velocity and velocity dispersion in source plane for individual bins. Velocity was extracted along the major kinematic axis with an aperture of 0.5 arcsec both for the data (back circles) and for the models (thick coloured lines) and is inclination corrected. Velocity dispersion corresponds to the intrinsic velocity dispersion (beam-smearing corrected) and was extracted in the full galaxy. A vertical bar marks the effective radius (R_e).

Open in new tabDownload slide

Figure 7.

Resolved extinction and SFR density maps for AS1063-arc (top) and A370 (bottom). In each row, from left to right: attenuation (τ) map, SFR surface density, from H β (units of [M_⊙ yr⁻¹ kpc⁻²]) and 1D radial profiles of both. These were obtained by averaging the values of extinction and SFR within distorted elliptical annuli of increasing radius (see Fig. 1). These distorted elliptical annuli correspond to circular apertures in source plane. from the centre of the galaxy. The position of the strongest star-forming clumps are marked with a black arrow.

Open in new tabDownload slide

Figure 8.

Velocity dispersion versus SFR density for several sub-galactic samples: in purple and blue circles, the results from AS1063 and A370, per Voronoi bin; in brown squares, the Genzel et al. (2011) giant star-forming clumps at z ∼ 2; in green triangles the sub-kiloparsec clumps at 1 < z < 4 of Livermore et al. (2015) and in dark triangles the clumps from Swinbank et al. (2012b). The lines show the results of the fit. Purple and blue full lines are fits to AS1063-arc and A370-sys1, respectively, dotted line the fit only to the Livermore et al. (2015), Swinbank et al. (2012b), and Genzel et al. (2011), and dashed dotted line fit to all data.

Open in new tabDownload slide

Figure 9.

[O ii ] surface brightness and velocity dispersion for the full sample except M0416-sys28 where not enough data was available. The Spearman coefficients of these two quantities are plot in the upper left corner of each panel. Unreliable points, (in black circles in the plots, see the text for details) were not included in the calculation of the Spearman coefficients. There is no correlation between [O ii ] surface brightness and velocity dispersion in the sample.

Open in new tabDownload slide

Finally, here we only focus on the ionized gas, traced by emission lines. The analysis of the properties of molecular gas, made from high-resolution ALMA data obtained for a number of these objects, will allow to probe the turbulence of the cold gas and to have a more complete view of the process of formation of star-forming clumps.

ACKNOWLEDGEMENTS

We warmly thank Benjamin Johnson for helping using prospector and updating the code to fit spectra and photometry simultaneously. We thank Jorge Sánchez Almeida and Andrea Ferrara for insightful suggestions about clumps, turbulence, and energy. We thank Mark Swinbank, Philippe Amram, Miroslava Dessauges, and Mathieu Puech for useful comments on this work. We thank David Bina for reducing A2390 and A2667 MUSE data. We also thank the referee for a constructive report with useful suggestions.

VP, DC, JR, JM, DL, and BC acknowledge support from the ERC starting grant 336736-CALENDS. TC, NB, and BE acknowledge support from the ANR FOGHAR (ANR-13-BS05-0010-02) and the OCEVU Labex (ANR-11-LABX-0060). TC acknowledges support from the A*MIDEX project (ANR-11-IDEX-0001-02) funded by the ‘Investissements d'avenir’ French government programme. KBS was supported by the BMBF grant 05A14BAC. JB acknowledges support by Fundação para a Ciência e Tecnologia (FCT) through national funds (UID/FIS/04434/2013), by FEDER through COMPETE2020 (POCI-01-0145-FEDER-007672), and by FCT through Investigador FCT contract IF/01654/2014/CP1215/CT0003. JR, FB, and RP also acknowledge support by the Programa de Cooperación Científica ECOS SUD Program C16U02. FB acknowledges support by CONICYT grants Basal-CATA PFB-06/2007, FONDECYT Regular 1141218, and the Ministry of Economy, Development, and Tourism's Millennium Science Initiative through grant IC120009, awarded to The Millennium Institute of Astrophysics, MAS.

Based on observations made with ESO Telescopes at the La Silla Paranal Observatory under programme IDs 060.A-9345, 094.A-0115, 095.A-0181, 096.A-0710, 097.A-0269, 100.A-0249, and 294.A-5032. Also based on observations obtained with the NASA/ESA Hubble Space Telescope, retrieved from the Mikulski Archive for Space Telescopes (MAST) at the Space Telescope Science Institute (STScI). STScI is operated by the Association of Universities for Research in Astronomy, Inc. under NASA contract NAS 5-26555. Part of this work was supported by the ECOS SUD Program C16U02. This research made use of astropy, a community-developed core python package for Astronomy (Astropy Collaboration 2013) and of VizieR catalogue access tool, CDS, Strasbourg, France.

Footnotes

REFERENCES