11institutetext:  Univ Lyon, Univ Lyon1, Ens de Lyon, CNRS, Centre de
Recherche Astrophysique de Lyon UMR5574, F-69230, Saint-Genis-Laval, France
22institutetext: Instituto de Astrofísica e Ciências do Espaço, Universidade
do Porto, CAUP, Rua das Estrelas, PT4150-762 Porto, Portugal 33institutetext:
Aix Marseille Université, CNRS, CNES, LAM (Laboratoire d’Astrophysique de
Marseille) UMR 7326, 13388, Marseille, France 44institutetext: Department of
Astronomy, University of Wisconsin, 475 N. Charter Street, Madison, WI 53706,
USA 55institutetext: Centre for Astrophysics and Supercomputing, Swinburne
University of Technology, Melbourne, VIC 3122, Australia 66institutetext:
Leibniz-Institut für Astrophysik Potsdam (AIP), An der Sternwarte 16, 14482
Potsdam, Germany 77institutetext: Institut für Physik und Astronomie,
Universität Potsdam,Karl-Liebknecht-Str. 24/25, D-14476 Golm, Germany
88institutetext: Univ. Grenoble Alpes, CNRS, Grenoble INP, GIPSA-lab, 11 rue
des Mathématiques, Grenoble Campus BP46, F-38402, Saint Martin D’heres Cedex,
France 99institutetext: Laboratoire Lagrange, CNRS, Université Côte d’Azur,
Observatoire de la Côte d’Azur, CS 34229, 06304, Nice, France
1010institutetext: Department of Astronomy, University of Texas at Austin,
2515 Speedway, Austin, TX 78712, USA 1111institutetext: INAF – Osservatorio di
Astrofisica e Scienza dello Spazio di Bologna, Via P. Gobetti 93/3, 40129
Bologna, Italy 1212institutetext: Observatoire de Genève, Université de
Genève, 51 Ch. des Maillettes, CH-1290 Versoix, Switzerland 1313institutetext:
Leiden Observatory, Leiden University, P.O. Box 9513, 2300 RA Leiden, The
Netherlands 1414institutetext: IRAP, Institut de Recherche en Astrophysique et
Planétologie, CNRS, Université de Toulouse, 14, avenue Edouard Belin, F-31400
Toulouse, France 1515institutetext: European Southern Observatory, Av. Alonso
de Córdova 3107, 763 0355 Vitacura, Santiago, Chile 1616institutetext: Max
Planck Institute for Astronomy, Königstuhl 17, 69117, Heidelberg, Germany
1717institutetext: Institut für Astrophysik, Universität Göttingen, Friedrich-
Hund-Platz 1, D-37077 Göttingen, Germany 1818institutetext: Canada-France-
Hawaii Telescope, CNRS, 96743 Kamuela, Hawaii, USA

# The MUSE Hubble Ultra Deep Field surveys: Data release II ††thanks: Based on
observations made with ESO telescopes at the La Silla Paranal Observatory
under the programs 094.A-0289(B), 095.A-0010(A), 096.A-0045(A), 096.A-0045(B)
and 1101.A-0127

Roland Bacon 11 Jarle Brinchmann 22 Simon Conseil 33 Michael Maseda 44 Themiya
Nanayakkara 55 Martin Wendt 6677 Raphael Bacher 1188 David Mary 99 Peter M.
Weilbacher 66 Davor Krajnović 66 Leindert Boogaard 1616 Nicolas Bouché 11
Thierry Contini 1414 Benoit Epinat 331818 Anna Feltre 1111 Yucheng Guo 11
Christian Herenz 1515 Wolfram Kollatschny 1717 Haruka Kusabe 1212 Floriane
Leclercq 1010 Léo Michel-Dansac 11 Roser Pello 33 Johan Richard 11 Martin Roth
66 Gregory Salvignol 11 Joop Schaye 1313 Matthias Steinmetz 66 Laurence Tresse
33 Tanya Urrutia 66 Anne Verhamme 1212 Eloise Vitte 1212 Lutz Wisotzki 66
Sebastiaan L. Zoutendijk 1313

(2022-10-15 Accepted)

We present the second data release of the MUSE Hubble Ultra-Deep Field
surveys, which includes the deepest spectroscopic survey ever performed. The
MUSE data, with their 3D content, amazing depth, wide spectral range, and
excellent spatial and medium spectral resolution, are rich in information.
Their location in the Hubble ultra-deep field area, which benefits from an
exquisite collection of ancillary panchromatic information, is a major asset.
This update of the first release incorporates a new 141-hour adaptive-optics-
assisted MUSE eXtremely Deep Field (MXDF; 1 arcmin diameter field of view) in
addition to the reprocessed 10-hour mosaic ($\rm 3\times 3\,arcmin^{2}$) and
the single 31-hour deep field ($\rm 1\times 1\,arcmin^{2}$). All three data
sets were processed and analyzed homogeneously using advanced data reduction
and analysis methods. The $3\sigma$ point-source flux limit of an unresolved
emission line reaches $\rm 3.1\times 10^{-19}$ and $\rm 6.3\times
10^{-20}\,\mathrm{erg\,s^{-1}\,cm^{-2}}$ at 10- and 141-hour depths,
respectively. We have securely identified and measured the redshift of 2221
sources, an increase of 41% compared to the first release. With the exception
of eight stars, the collected sample consists of 25 nearby galaxies
($z<0.25$), 677 [O ii] emitters ($z=0.25-1.5$), 201 galaxies in the MUSE
redshift desert range ($z=1.5-2.8$), and 1308 $\mathrm{Ly}\alpha$ emitters
($z=2.8-6.7$). This represents an order of magnitude more redshifts than the
collection of all spectroscopic redshifts obtained before MUSE in the Hubble
ultra-deep field area (i.e., 2221 versus 292). At high redshift ($z>3$), the
difference is even more striking, with a factor of 65 increase (1308 versus
20). We compared the measured redshifts against three published photometric
redshift catalogs and find the photo-z accuracy to be lower than the
constraints provided by photo-z fitting codes. Eighty percent of the galaxies
in our final catalog have an HST counterpart. These galaxies are on average
faint, with a median AB F775W magnitude of 25.7 and 28.7 for the [O ii] and Ly
emitters, respectively. Fits of their spectral energy distribution show that
these galaxies tend to be low-mass star-forming galaxies, with a median
stellar mass of $\rm 6.2\times 10^{8}\,M_{\odot}$ and a median star-formation
rate of $\rm 0.4\,M_{\odot}yr^{-1}$. We measured the completeness of our
catalog with respect to HST and found that, in the deepest 141-hour area, 50%
completeness is achieved for an AB magnitude of 27.6 and 28.7 (F775W) at
$z=0.8-1.6$ and $z=3.2-4.5$, respectively. Twenty percent of our catalog, or
424 galaxies, have no HST counterpart. The vast majority of these new sources
are high equivalent-width $z>2.8$ $\mathrm{Ly}\alpha$ emitters that are
detected by MUSE thanks to their bright and asymmetric broad Ly line. We
release advanced data products, specific software, and a web interface to
select and download data sets.

###### Key Words.:

Galaxies: high-redshift – Galaxies: evolution – Cosmology: observations –
Techniques: integral field spectroscopy

##  1 Introduction

In 2017 we published the first MUSE spectroscopic survey data release (Bacon
et al., 2017; Inami et al., 2017, hereafter DR1) in the Hubble ultra-deep
field (Beckwith et al., 2006) area. These two papers were published together
with eight other papers, providing a first glimpse into the richness of the
data and their scientific impact. In that series of papers we investigated the
photometric redshift properties of the sample (Brinchmann et al., 2017), the
properties of C iii] emitters (Maseda et al., 2017), the spatially resolved
stellar kinematics of galaxies at z=0.2–0.8 (Guérou et al., 2017), the faint
end of the Ly luminosity function (LF; Drake et al., 2017), the properties of
Fe ii∗ emission (Finley et al., 2017), the extended Ly halos (Leclercq et al.,
2017), the evolution of the galaxy merger fraction (Ventou et al., 2017), and
the equivalent-width properties of $\mathrm{Ly}\alpha$ emitters (Hashimoto et
al., 2017).

These initial studies did not exhaust the scientific content of the data set
and were quickly followed by many others: the study of Mg ii emission and
absorption in star-forming galaxies (Feltre et al., 2018), the recovery of
systemic redshifts from Ly line profiles (Verhamme et al., 2018), the covering
fraction of Ly emission (Wisotzki et al., 2018), the low-mass end of the star
formation sequence (Boogaard et al., 2018), the spatially resolved properties
of Ly halos (Leclercq et al., 2020), the evolution of the $\mathrm{Ly}\alpha$
emitter fraction (Kusakabe et al., 2020), the molecular gas properties of high
redshift galaxies (Boogaard et al., 2019, 2020, 2021; Inami et al., 2020), the
study of high equivalent-width $\mathrm{Ly}\alpha$ emitters (Maseda et al.,
2018, 2020), the study of He ii emission line properties (Nanayakkara et al.,
2019), the rest-UV properties of $\mathrm{Ly}\alpha$ emitters (Feltre et al.,
2020), and the angular momentum of low-mass star-forming galaxies (Bouché et
al., 2021).

Shortly after the completion of the MUSE Hubble Ultra-Deep Field (hereafter
HUDF) survey, MUSE was successfully coupled to the ground-layer module of the
VLT adaptive optics facility (AOF; Leibundgut et al., 2017; Madec et al.,
2018). The new system offers improved spatial resolution as well as faster
survey times by allowing a wider range of atmospheric conditions (i.e., seeing
and airmass) for operations.

The DR1 was based on two data sets: a $\rm 3\times 3\,arcmin^{2}$ mosaic of
nine MUSE fields at a 10-hour depth (hereafter MOSAIC) and a single $\rm
1\times 1\,arcmin^{2}$ 31-hour depth field (hereafter UDF-10). The improved
depth of UDF-10 proved to be essential, especially for the study of the faint
end of the $\mathrm{Ly}\alpha$ emitter population and the study of diffuse
emission. We then decided to push forward in depth and start a new adaptive-
optics-assisted GTO survey in the same area with the goal of reaching a depth
of over 100 hours. This new survey, called the MUSE eXtremly Deep Field
(MXDF), was completed in January 2019. With an achieved depth of 141 hours, it
is a key addition to the existing MOSAIC and UDF-10 data sets and is the
deepest spectroscopic survey ever performed.

In addition to these new observations, we significantly improved the data-
reduction and data-analysis system with respect to DR1 (Sect. 6.1.2). The
MOSAIC and UDF-10 data sets were therefore reprocessed with the same tools and
methodology used for the MXDF processing. This second data release (DR2),
which incorporates these three data sets, aims to provide a comprehensive and
homogeneous deep spectroscopic survey in the HUDF iconic field.

The paper is organized as follows: the MXDF observations are reported in Sect.
2 and the improved data reduction in Sect. 3 and Appendix B. Section 4
presents the data properties. The processes of source detection and
classification are described in Sect. 5, and the resulting catalogs and sample
properties are presented in Sect. 6. A summary and conclusions are given in
Sect. 7. Finally, the released advanced data products and software are
detailed in Appendix A.

Data reduction and catalog building for deep integral field unit (IFU) data is
not yet common. Therefore, we aim to provide sufficient details on the
methodology used in all our data processing and data analysis steps, to enable
future users to perform similar deep fields. Readers not interested in these
technical details can skip Sects. 3, 4, and 5.

##  2 Observations

In this section we report only recent observations of MXDF (GTO Large Program
1101.A-0127, PI R. Bacon). For previous observations related to the UDF-10 and
MOSAIC fields, we refer to Bacon et al. (2017). The observing campaign started
in August 2018 and lasted until January 2019 for a total of 6 runs made during
new moon periods.

All observations were made with the VLT’s AOF and GALACSI, its dedicated
ground-layer adaptive optics (GLAO) system (Kolb et al. 2016; Madec et al.
2018). With respect to non-AOF observations, the only change in the MUSE
instrumental configuration is the notch filter that blocks bright light due to
the four sodium laser guide stars in the 5800-5966 Å wavelength range (Fig.
1). The AOF ran smoothly and achieved robust performance during all runs.

The location of the MXDF field (Fig. 2) was chosen (i) to be in the HUDF
extremely deep field region with WFC3 deep imaging (Illingworth et al., 2013)
(ii) to match the AOF tip/tilt natural star requirement (i.e., V $>$ 18.5 in a
3.5 arcmin field of view and outside the MUSE field of view), and (iii) to
have a usable slow-guiding star in the outer circle to compensate for the
derotator’s wobble and possible misalignment of the Nasmyth platform with
respect to the telescope’s focal plane.

As shown in Table 1, we observed the field over a wide range of atmospheric
conditions, including poor seeing (1.2 arcsec) and up to high air mass (1.4).
The ground layer fraction, a key parameter for GLAO observations – the larger
the better – averaged 65%, and rarely was less than 50%. We observed some
correlation between seeing and ground layer fraction: on many occasions, poor
seeing is related to increased ground layer turbulence. This behavior of the
Paranal atmospheric turbulence is fortunate because it makes the GLAO mode
very effective.

It is instructive to compare the observing conditions of the UDF-10 and MOSAIC
campaigns with the MXDF campaign (Table 1). The first campaign, conducted
before the AOF was commissioned, was limited to the natural seeing observing
mode of MUSE. To maximize the final spatial resolution, we therefore decided
to observe the UDF field only under good seeing conditions (seeing less than
0.8 arcsec) and at low air mass (less than 1.2). The consequence is that it
took 2.5 years and 8 runs to accumulate the required telescope time. In
contrast, the MXDF campaign took only 6 months and 6 consecutive runs to reach
155 hours of integration time. As we will see in Sect. 4.2.1 (Fig. 5), thanks
to the performance of the AOF/GALACSI GLAO, the spatial resolution achieved by
MXDF is better than that of UDF-10 and MOSAIC, despite its relaxed atmospheric
observing conditions.

Table 1: Summary of MXDF, UDF-10, and MOSAIC observations.

Field | Runsaa$a$Starting and end dates. The number of conducted runs is indicated in parentheses. | Integbb$b$Total open shutter time in hours. | Airmasscc$c$50% (q50, median), 10% (q10) and 90% (q90) percentile. | Seeingc,dc,dfootnotemark: $c,d$ (″) | GL fracc,ec,efootnotemark: $c,e$
---|---|---|---|---|---
|  |  | q50 | q10 | q90 | q50 | q10 | q90 | q50 | q10 | q90
MXDF | 08/18-01/19 (6) | 155 | 1.08 | 1.01 | 1.37 | 0.76 | 0.54 | 1.11 | 0.65 | 0.46 | 0.83
|  UDF-10 &
---
MOSAIC
09/14-02/16 (8) | 116 | 1.05 | 1.00 | 1.19 | 0.75 | 0.60 | 0.86 |  |  |

111

The MXDF observing strategy also differs from the 90-degree rotation plus
small offset dithering patterns commonly used in MUSE observations. Each
observing block is similar to the scheme used in the first campaign: that is,
a set of four 25 min exposures with successive 90∘ instrument rotations. But
in order to further reduce the systematics and to break the horizontal and
vertical patterns introduced by the instrument field splitter and slicer
geometry, we systematically rotated the field of view by a few degrees between
each observing block. Consequently, the final combined field of view is
approximately circular (Fig. 2) with a radius of 41″ and 31″ for respectively
10+ and 100+ hours of depth. The field center celestial coordinates are
53.16467 deg., -27.78537 deg. (J2000 FK5). After the rejection of a few bad
exposures due to satellite track contamination, poor final spatial resolution,
or cloud absorption, the achieved final maximum depth was 141 hours.

Figure 1: MUSE spectral range. The typical sky spectrum (shown in black) is
observed in one 25 min observing block. The blue shaded wavelength region
shows the location of the sodium notch filter used in MUSE GLAO mode. The
response curves of the HST ACS filters F606W, F775W, and F850LP are also
indicated.  Figure 2: Location of the three deep fields used in this paper:
MXDF (141-hour depth), MOSAIC (10-hour depth), and UDF-10 (31-hour depth)
overlaid on the HST F775W UDF image. The dotted and dashed red circles show
the MXDF 10- and 100-hour exposure time isocontours, respectively.

##  3 Data reduction

After DR1, described in Bacon et al. (2017), we continued to work on the data
reduction process, which led to several major improvements. The overall
process is similar to DR1, with important changes in the self-calibration
algorithm, in the sky-subtraction with the Zurich Atmospheric Purge (ZAP)
software (Soto et al. 2016, Appendix B.3), and the use of a “superflat.” This
process was first applied to the UDF-10 and MOSAIC (278 exposures) and then to
the MXDF (373 exposures). The resulting depth and color images for the three
data sets are shown in Fig. 3 and Fig. 4, respectively. A detailed description
of this improved data reduction process is given in Appendix B.

Figure 3: Achieved depth in hours. The MOSAIC, UDF-10, and MXDF exposure maps
are shown in the first three panels, starting from the left. The Max exposure
map, computed as the maximum depth for each spaxel, is shown in the right
panel. UDF-10 and MXDF are to the north (top) and the east (left),
respectively, while the MOSAIC and Max exposure maps are rotated by 42∘(see
Fig. 2). Axis labels are in arcmin.  Figure 4: Reconstructed pseudo-color
images of the MOSAIC, MXDF, and UDF-10 data sets.

##  4 Data properties

###  4.1 Astrometry and photometry

As in DR1, the datacube’s world coordinate system has been matched to the
Hubble ACS astrometry (Appendix B). By construction, the astrometry of the
datacube should therefore be consistent with the HUDF published photometric
catalogs: UVUDF (Rafelski et al., 2015), 3D-HST (Skelton et al., 2014), and
CANDELS v2 (Whitaker et al., 2019). We note that Dunlop et al. (2017) and
Franco et al. (2018) have shown that the Hubble Space Telescope (HST)
astrometry is offset from the Gaia DR2 catalogs (Gaia Collaboration et al.,
2016, 2018). The offset, as measured by Whitaker et al. (2019), amounts to
$\Delta\text{Ra}=+0.094\pm 0.042$ arcsec and $\Delta\text{Dec}=-0.26\pm 0.10$
arcsec. In this paper, we have ignored this offset in order to maintain the
same astrometric reference as the HST catalogs.

We measured the astrometric accuracy achieved by comparing source centroids in
the HST ACS F775W high resolution image with the corresponding values measured
in MUSE reconstructed image. We proceeded as follows: For each field, we
derived the equivalent F775W broadband image by computing at each spaxel the
transmission-weighted average flux of the corresponding F775W ACS filter (Fig.
1). This reconstructed image is then convolved by the HST ACS point spread
function (PSF). We approximated this PSF using a Moffat function222The
circular 2D Moffat profile (Moffat, 1969) has two parameters: the FWHM in
arcsec and the $\rm\beta$ shape parameter, $\rm\beta=1$ and $\rm\beta\gg 1$
for respectively a Lorentzian or a Gaussian distribution. with a 0.085 arcsec
full width at half maximum (FWHM) and $\beta=1.60$. The HST ACS F775W image is
convolved with the MUSE PSF value at the F775W reference wavelength (7750 Å),
that is, 0.45 arcsec FWHM and $\beta=1.89$, 0.60 arcsec FWHM and $\beta=2.80,$
and 0.63 arcsec FWHM and $\beta=2.80$ for the MXDF, UDF-10, and MOSAIC data
sets, respectively. Source detection is performed on both convolved images
with SExtractor (Bertin & Arnouts, 1996). We tune SExtractor’s detection
parameters to obtain similar segmentation maps between HST and MUSE convolved
images. The two catalogs are crossmatched and the mean offset and its standard
deviation are derived using $3\sigma$ sigma-clip statistics on the matched
catalog. The astrometric differences between the MUSE datacubes and HST ACS
F775W image are given in the upper part of Table 2. If we exclude the faint
sources with $\text{AB}>27$, the mean astrometric relative error is about 0.1
arcsec, which is half of a spaxel or about 1/6 of the spatial resolution. We
note that astrometric errors are larger for MXDF than for the MOSAIC data set.
Such an increase with depth is expected given the increase of source
confusion.

We used the same crossmatched catalogs to evaluate the achieved accuracy in
source AB magnitude for MUSE data sets with respect to HST photometry. The
results, reported in the lower part of Table 2, are based on SExtractor
computed automatic aperture magnitude333SExtractor’s automatic aperture
photometry routine is derived from Kron’s algorithm (Kron, 1980); see
https://sextractor.readthedocs.io/en/latest/Photom.html for more detailed
information.. We note that the magnitude difference is below 0.1 for MXDF
sources with $\text{AB}<27$. As expected, the magnitude differences decrease
with the source brightness and increase with depth.

Table 2: Estimated astrometric and photometric accuracy for the three data
sets when compared to the HST ACS F775W image. Mean and standard deviation are
given for three AB magnitude intervals.

data set | $\text{AB}<25$ | $25<\text{AB}<27$ | $27<\text{AB}$
---|---|---|---
| Average astrometric offsets (arcsec)
MXDF | $0.12\pm 0.04$ | $0.13\pm 0.05$ | $0.13\pm 0.06$
UDF10 | $0.08\pm 0.02$ | $0.10\pm 0.05$ | $0.15\pm 0.09$
MOSAIC | $0.05\pm 0.02$ | $0.08\pm 0.05$ | $0.16\pm 0.10$
| Average AB magnitude differences
MXDF | $0.08\pm 0.04$ | $0.09\pm 0.07$ | $0.25\pm 0.19$
UDF10 | $0.07\pm 0.04$ | $0.23\pm 0.13$ | $0.35\pm 0.24$
MOSAIC | $0.10\pm 0.08$ | $0.24\pm 0.16$ | $0.32\pm 0.23$
444The mean AB magnitude error reported by SExtractor for the HST convolved
image amounts to 0.002, 0.011 and 0.030 in the $\text{AB}<25$,
$25<\text{AB}<27$ and $\text{AB}>27$ intervals, respectively.

###  4.2 Spatial and spectral resolution

####  4.2.1 MXDF spatial PSF

Spatial resolution is an important ingredient for most analyses. The spatial
PSF was first estimated for each individual raw exposure to reject bad
exposures before the final datacube was produced. Then, the final PSF was
estimated on the combined datacube.

Unfortunately, there is no point source bright enough in the MXDF field of
view to estimate the PSF in a single exposure of 25 min integration. Moreover,
the seeing value reported by the Paranal seeing monitor is not sufficient by
itself to estimate the PSF given the adaptive optics (AO) correction provided
by the AOF GLAO system. We therefore used the tool muse-psfr developed by
Fusco et al. (2020) to derive an estimate of the PSF using AO telemetry.
During an exposure, atmospheric parameters such as seeing, ground-layer
turbulence fraction, turbulence outer scale, are derived from wavefront sensor
telemetry and regularly recorded in a table in the raw FITS exposure. The
algorithm uses this information to produce an estimate of the spatial PSF in
the form of a circular 2D Moffat profile. Fusco et al. (2020) measure a
standard deviation of 0.06 arcsec for the derived FWHM. We note that some
systematic offsets were found at high spatial resolution (Fig. 17 of Fusco et
al., 2020). We also found a similar offset when we compared the derived value
of muse-psfr with a direct fit to the combined datacube (Table 4).

We show in Fig. 5 the distribution of the computed MXDF FWHMs555A value of
0.06 arcsec was added to the values of MXDF muse-psfr to account for the
measured offset between muse-psfr and the direct Moffat fit. at 7000 Å for the
373 individual exposures. The benefit of the GLAO AOF system is evident when
the histogram of MXDF FWHMs is compared to the corresponding MOSAIC plus
UDF-10 values. Due to the specific rotation scheme used for the MXDF
observations, the outer ring of the field results from the combination of
different exposures and may thus have a potentially different PSF than the
central part. We show in Fig. 6 that this is indeed the case with a variation
of 0.10 arcsec FWHM. We note, however, that within the area with depth $>100$
hours there is no detectable spatial variation in the PSF.

Although there is no bright star in the field of view, there is a fainter M
star (ID 5102 in R15 catalog, F775W AB magnitude 24.7) located in the 100-hour
depth area that can be used to measure the PSF. To take into account the
evolution of the PSF with wavelength, we divided the final combined datacube
into 20 wavelength slices of 232 Å each and fit a circular Moffat function to
each image. As shown in Fig. 7 (left panels), the circular Moffat function
represents the data well. The figure also shows the evolution of the shape
parameters FWHM and $\rm\beta$ as a function of wavelength. Except at the blue
end of the MUSE wavelength range, where the M star becomes too faint, the
shape parameter $\rm\beta$ is approximately constant. We adopt
$\rm\beta=2.123$ and perform the Moffat fit again, leaving the FWHM as a
variable. The adopted third-order polynomial fitted to the evolution of the
FWHM as a function of wavelength is shown in the upper right panel of Fig. 7.
The final spatial PSF is then given by

$\displaystyle\rm PSF(r,\lambda)$
$\displaystyle=(\beta-1)/(\pi\alpha(\lambda)^{2})\times(1+r^{2}/\alpha^{2})^{-\beta}$
(1) $\displaystyle\rm\alpha(\lambda)$
$\displaystyle=FWHM(\lambda)\times\left[2\sqrt{2^{1/\beta}-1}\right]^{-1}$
$\displaystyle\rm FWHM(\lambda)$ $\displaystyle=a_{3}\lambda^{\prime
3}+a_{2}\lambda^{\prime 2}+a_{1}\lambda^{\prime}+a_{0}$
$\displaystyle\rm\lambda^{\prime}$
$\displaystyle=(\lambda-4850)/(9350-4850)-0.5.$

The wavelength ($\lambda$) and FWHM units are Å and arcsec, respectively. The
$\beta$ value and the FWHM polynomial coefficients ($a_{n}$) are given in
Table 3.

We compared this PSF estimate with the values derived by convolution of the
HST images for the F606W, F775W, F814W, and F850LP broadband filters. The
method is identical to that used in the paper I (Bacon et al., 2017, Sect.
5.1) for the MOSAIC and UDF-10 datacubes. For ease of comparison, a direct
Moffat fit was performed on the M star for the four reconstructed HST
broadband MUSE images. An additional comparison was performed by fitting a
Moffat distribution to the average of all individual muse-psfr PSFs. The
Moffat FWHM and $\beta$ shape parameter for the three different estimates are
shown in the Table 4. Since the FWHM and $\beta$ parameters are not
independent, we also provide the FWHM obtained by fitting a Gaussian profile
to judge the difference in FWHM. However, we recall that the MUSE PSF is not
well represented by a Gaussian function (Husser et al., 2016; Weilbacher et
al., 2020). As shown in Table 4, all methods give similar results with the
expected decrease in FWHM with increasing wavelengths. However, there is a
systematic offset between the methods: the HST convolution method
overestimates the FWHM by $\approx$+0.05 arcsec while the muse-psfr method
underestimates it by $\approx$-0.06 arcsec.

Figure 5: Spatial resolution of individual MXDF observations compared to the
UDF-10 plus MOSAIC observations. The FWHMs (in arcsec) are derived from the
Moffat PSF model at 7000 Å. The normalized FWHM histograms of the MXDF and
UDF-10 plus MOSAIC observations are shown in green and gray, respectively. The
FWHM MXDF values have been corrected for a $\rm+0.06$ arcsec offset (see
text). Figure 6: Field evolution of the MXDF spatial resolution. The colors
indicate the mean value (in arcsec) of all individual Moffat FWHMs at 7000 Å
at each spaxel location as measured by the PSF AO reconstruction algorithm and
corrected for an offset of $\rm+0.060$ arcsec. The red circle displays the
100-hour depth contour. Figure 7: MXDF Moffat fit of the M star RID-5102. The
fit is performed on 20 narrowband images of 232 Å width, equally distributed
along the MXDF cube wavelength axis. Examples of the first, central, and last
images (each box side is 5 arcsec in length) are displayed in the first
column. The corresponding Moffat fit and residuals are shown in the two other
columns. The FWHM and $\beta$ Moffat fitted values are shown in the upper and
bottom right panels, together with their polynomial approximation.  Table 3:
Circular Moffat PSF models for the three MUSE data sets.

data set | $\beta$ | $a_{0}$ | $a_{1}$ | $a_{2}$ | $a_{3}$ | $F_{4850}$ | $F_{7000}$ | $F_{9350}$
---|---|---|---|---|---|---|---|---
MXDF | 2.123 | 0.5465 | -0.2474 | 0.1975 | -0.0295 | 0.72 | 0.55 | 0.47
UDF-10 | 2.800 | 0.6179 | -0.1353 |  |  | 0.69 | 0.62 | 0.55
MOSAIC | 2.800 | 0.6378 | -0.1515 |  |  | 0.71 | 0.64 | 0.56
666 $\beta$ is the Moffat shape parameter and $a_{n}$ are the coefficients of
the polynomial approximation to FWHM($\lambda$) defined in Eq. 1. In the case
of the MOSAIC data set, we provide the median PSF for the nine fields. FWHM
values in arcsec (e.g., $F_{4850}$) are given for the blue (4850 Å), central
(7000 Å) and red (9350 Å) wavelengths.

Table 4: Comparison of different spatial MXDF PSF estimates for four HST
reconstructed broadband images.

Filter | Moffat FWHM | Moffat $\rm\beta$ | Eq. Gaussian FWHM
---|---|---|---
| Star | HST | AO tel | Star | HST | AO tel | Star | HST | AO tel
F606W | 0.60 | 0.64 | 0.54 | 2.10 | 2.41 | 2.01 | 0.72 | 0.75 | 0.65
F775W | 0.52 | 0.56 | 0.45 | 2.17 | 2.25 | 1.88 | 0.61 | 0.66 | 0.55
F814W | 0.50 | 0.55 | 0.43 | 2.19 | 2.35 | 1.84 | 0.59 | 0.64 | 0.53
F850LP | 0.48 | 0.51 | 0.41 | 2.25 | 2.25 | 1.80 | 0.57 | 0.60 | 0.50
777 FWHM (arcsec) and $\rm\beta$ values are given for the direct Moffat fit of
the RID-5102 M star (Star columns), the comparison with HST (HST columns) and
the mean of the PSF reconstruction values (AO tel columns). The Gaussian
equivalent FWHMs are given in the last three columns.

####  4.2.2 Line spread function

As in DR1, we measure the line spread function (LSF) on a combined cube of
data where the sky has not been subtracted. A Gaussian fit was performed on
the bright sky lines taking into account the contribution of adjacent fainter
lines. As shown in Fig. 8, the MXDF LSF shows the expected instrumental
evolution with wavelength. A good model of the median LSF FWHM is given by the
UDF-10 model: $\rm LSF(\lambda)=5.866\times 10^{-8}\lambda^{2}-9.187\times
10^{-4}\lambda+6.040$, with LSF and $\lambda$ in Å. The LSF is approximately
constant in the field of view, except in the outer rings, which show larger
variation. This is mainly due to the limited number of exposures that are
combined together at the edges, which leaves the intrinsic variation between
the 24 IFUs more visible.

Figure 8: Measurements of the MXDF LSF in the combined datacube (not sky
subtracted). The bottom panel displays the fitted Gaussian FWHM for a few
bright sky lines. Variations with the field of view are given by the error
bars (standard deviation). The DR1 UDF-10 polynomial approximation is shown as
a dashed line. Examples of FWHM spatial maps are given for three sky lines in
the top panels. FWHM units is Å.

###  4.3 Noise property and limiting flux

To compare the noise level of the three data sets, we selected 214 one-
arcsecond-diameter apertures located in empty locations and common to all data
sets. The location of the aperture was chosen to have a similar depth for each
data set: $140.8\pm 0.06$ hours, $31.2\pm 0.2$ hours, and $10.0\pm 0.3$ hours
for the MXDF, UDF-10, and MOSAIC data sets, respectively. The corresponding
noise level as a function of wavelength is calculated as the median of the 214
standard deviations. It should scale with exposure time ($t$) as
$\rm\sqrt{t}$. In Fig. 9 we display the UDF-10 and MXDF noise level scaled by
the square root of the depth and divided by the MOSAIC scaled value. The
figure shows that the deviation from the theoretical expectation is indeed
very small (i.e., 10-20%). The slight improvement with depth results from a
better control of the systematics due to the observation strategy (dithering
and rotation; see Sect. 2).

We estimated the 1$\sigma$ surface brightness limit for an unresolved lines
($\textit{SB}_{line}$) with the following formula:

$\textit{SB}_{line}(\lambda)=\frac{\sqrt{V_{line}(\lambda)}\,\Delta\lambda}{\sqrt{25}\,\Delta
s^{2}}\\\ \mathrm{with}\\\
V_{line}(\lambda)=\sum_{\lambda^{\prime}=\lambda_{0}-k\sigma(\lambda)}^{\lambda_{0}+k\sigma(\lambda)}V_{s}(\lambda^{\prime}),$
(2)

where $V_{s}(\lambda)$ is the variance by spaxel, $\Delta\lambda$ the
wavelength bin size (1.25 Å), $\Delta s$ the spaxel bin size (0.2 arcsec), 25
the number of spaxels in 1 arcsec2 and $\sigma(\lambda)$ is derived from the
LSF approximation given in Sect. 4.2.2 with
$\sigma(\lambda)=\text{LSF}(\lambda)/2.355$. The wavelength interval
($k=1.29$) was chosen to capture 90% of the line flux. The variance by spaxel
was computed on the central MXDF area with a depth larger than 120 hours after
masking the bright sources and taking the median of the datacube variance. We
note that the variance is already corrected for the correlated noise as
described in Appendix B.6.

The same computation was performed for the UDF-10 data set, restricted to the
area with depths greater than 30 hours, and the MOSAIC data set. Results are
shown in Fig. 10. The best sensitivity is achieved in the MXDF deep area and
in the 6800-7800 Å wavelength range (outside sky lines) where an unresolved
emission line surface brightness of $\mathrm{1.0\times
10^{-19}\,erg\,s^{-1}\,cm^{-2}\,arcsec^{-2}}$ will be detected at 3$\sigma$.

We perform a similar estimation for a point-like source with an unresolved
emission line, using the spatial PSF estimated in Sect. 4.2.1 with the
following formula:

$F_{pl}(\lambda)=\Delta\lambda\,\sqrt{\frac{\pi r^{2}}{\Delta
s^{2}}V_{line}(\lambda)},$ (3)

where $V_{line}$ is defined in Eq. 2 and $r$ is the 80% enclosed flux circular
radius for the corresponding Moffat PSF model at 7000 Å. The best point-source
sensitivity (Fig. 11) is achieved in the MXDF central area and in the
6700-9300 Å wavelength range (sky lines excluded) with a 3$\sigma$ detection
limit of $6.3\times 10^{-20}\mathrm{erg\,s^{-1}\,cm^{-2}}$ for an unresolved
emission line.

Figure 9: Scaled raw noise comparison between the three data sets. The plot
displays the UDF-10 and MXDF to MOSAIC scaled noise ratio. The scaled noise is
computed as the median noise level for each data set in a common set of 1″
empty apertures, after scaling by the square root of the depth in hours. We
note that the noise level has not been corrected for correlations caused by
the interpolation.  Figure 10: Unresolved emission line median 1$\sigma$
surface brightness limit for the three data sets in units of
$\mathrm{erg\,s^{-1}\,cm^{-2}\,arcsec^{-2}}$.  Figure 11: Unresolved emission
line median 1$\sigma$ point-like source limit for the three data sets in units
of $\mathrm{erg\,s^{-1}\,cm^{-2}}$.

##  5 Source detection and classification

The next stage after the production of the final datacubes for the three data
sets (MOSAIC, UDF-10 and MXDF) is the detection and classification of sources.
This is a delicate task, involving many steps to ensure the data quality and
homogeneity of the final product. An overview, followed by detailed
descriptions of each of the step, is presented in the following subsections.

###  5.1 Overview

A schematic of the process is shown in Fig. 12. The first phase is the
production of inspection data sets (process chart one in Fig. 12). Candidate
sources are produced in two parallel branches: (i) from blind detection using
ORIGIN software (Sect. 5.2) and (ii) from de-blended extraction using the
ODHIN algorithm (Sect. 5.3) based on the R15 HST catalog and segmentation map.
For each candidate source, the first five redshift solutions obtained with the
pyMarZ software are recorded (Sect. 5.4). For each redshift solution, the
corresponding emission and absorption lines are computed with pyPlatefit
(Sect. 5.5) and narrow bands for the highest S/N lines are derived (Sect.
5.6). The data resulting from the different processes are organized into a
database with a set of candidate source catalogs, one for each input datacube
and each detection branch (i.e., ORIGIN and ODHIN).

The second phase (process chart two in Fig. 12, Sect. 5.7) consists of the
visual inspection of a set of candidate sources by a group of experts who
select the appropriate redshift solution, assigns confidence level and
crossmatch the source with the HST R15 and DR1 catalogs. Conflicts between the
experts’ solutions are identified and resolved in reconciliation meetings.
Duplicates between ORIGIN and ODHIN sources are sorted out and the retained
solutions are attributed a unique MUSE identifier888We use the following
naming convention: MID-nnnn for MUSE identifiers and RID-nnnn for R15
identifiers.

In the last step (process chart three in Fig. 12, Sect. 5.8), we selected the
most appropriate extraction and refined the redshift solution, emission, and
absorption line fits and narrow bands for each source. The results are
inspected and corrective actions are identified and executed. This process is
iterated until the results are deemed satisfactory. At the end of this last
step, the final data products are produced, that is, a set of catalogs
(Appendix A.2) and data organized in the MUSE Python Data Analysis Framework
(MPDAF) source format (Appendix A.3).

Figure 12: Schematic of the processes involved in the data release production.

###  5.2 ORIGIN blind source detection

The ORIGIN software (Mary et al., 2020) has been developed to automatically
detect faint line emitters in MUSE datacubes. It is optimized for the
detection of compact sources with faint spatial-spectral emission signatures
and provides an automated and reliable estimate of the purity (i.e., one minus
the fraction of false positives, Benjamini & Hochberg 1995; Mary & Roquain
2021). All three data sets were processed with the same input parameters.
Below we present the example of the MXDF data set.

The preprocessing steps (Fig. 13) perform a first pass continuum removal using
a discrete cosine transformation of order ten along the datacube wavelength
axis. The resulting white light continuum image is then segmented to mask the
bright continuum sources. A second segmentation is also performed on the S/N
residual image (i.e., the reconstructed continuum-free image of the datacube
divided by the square root of the variance). The two segmentation masks are
merged. The unmasked region defines the ”faint region” that will be used in
the following steps to avoid disturbance due to bright or noisy sources.

Although most of the continuum signal was removed in the preprocessing steps,
there are still residuals in the resulting datacubes that will create too many
false positives in the detection step, forcing the use of a higher detection
threshold and preventing the detection of faint emission line sources. The
continuum residuals, as well as the systematics present in the datacube, are
subsequently removed with the iterative Principal Component Analysis (PCA)
process described in Sect. 3.1.3 of Mary et al. (2020). The MXDF field of view
is segmented into six zones (Fig. 14 left panel; see Sect. 3.1.2 of Mary et
al. 2020). In each segment, a threshold is calculated on the spaxel
distribution to evaluate the deviation from the normal distribution999The
noise distribution in the MXDF is well represented by a normal distribution,
as are the UDF-10 and MOSAIC noise statistics; see Fig. 18 of Bacon et al.
2017.. This threshold is then used by the iterative PCA process to clean the
signal from nuisance sources. In the left panel of Fig. 15, we show the number
of cleaning iterations at each spaxel. We note that while in most of the
spaxels only a few iterations were needed to clean the data, many more were
required at specific locations (e.g., at the edge of the field or in the
presence of complex bright sources residuals).

The datacube is then filtered using a library of possible 3D signatures,
constructed as spectral line profiles spread spatially by the PSF and a
generalized likelihood ratio (GLR) is computed for each voxel, as described in
Sect. 3.2.1 of Mary et al. (2020). The image of the maxima of the resulting
GLR datacube is shown in Fig. 15, right panel. A comparison of the local
maxima counts in the GLR datacube, restricted to the faint area, with the
corresponding local minima counts is used to derive the purity as a function
of the threshold applied to the cube of local maxima (Fig. 16, left panel).
Using a purity of 0.8 for the GLR datacube, we identify 984 sources after
spatial grouping of the detected emission lines. A similar operation is
performed in the S/N datacube, prior to PCA subtraction, to recover possible
bright sources affected by the PCA process. In that case we use a higher
purity of 0.9 to search for these additional bright sources. This yields 18
more sources.

For each source, an optimal extraction is performed on the MUSE datacube with
the Horne (1986) algorithm and using as a weighting and segmentation image,
the image obtained by summing the GLR datacube over a wavelength window
centered on the detection peak wavelength and with a width of
$\pm\,2\times\text{FWHM}$, where FWHM is the width of the spectral template
that provides the highest correlation peak (see Sect. 3.5 of Mary et al.
2020). We emphasize that the extraction is always performed on the original
MUSE datacube, the ORIGIN pseudo narrowband image being used only as weighting
image.

The process was similar for the UDF-10 and MOSAIC data sets. We compared our
results with another blind source detection software program, LSDCat,
developed by Herenz & Wisotzki (2017) in the context of the MUSE-Wide survey
(Herenz et al., 2017; Urrutia et al., 2019; Schmidt et al., 2021). LSDCat is
based on a 3D matched filter, which is, to a first approximation, similar to
the mathematical basis of the ORIGIN algorithm. However, there are some
important differences (e.g., the continuum subtraction or the use of GLR) and
we should then be able to assess the robustness of our results by comparing
the two detection catalogs. We use the same MOSAIC datacube to perform this
comparison. Excluding the UDF-10 footprint to have a homogeneous depth, the
LSDCat catalog is composed of 1190 sources. This must be compared to the 1287
ORIGIN sources with a redshift confidence $\rm ZCONF>0$ (see Sect. 5.7.3)
using the same selection. 81 LSDCat sources (6%) are unmatched in the ORIGIN
catalog. A detailed examination of these sources shows that most of them are
low purity sources that were detected by ORIGIN but not confirmed in the
inspection process. Of the ORIGIN sources, 246 (19%) are not matched to the
LSDCat catalog. Half of these sources are low confidence ($\rm ZCONF=1$)
sources that fall below the LSDCat detection threshold. Overall, the vast
majority of high confidence sources (93%) are detected by the two methods,
giving confidence in the robustness of the results.

Figure 13: ORIGIN preprocessing steps. From left to right, the MXDF MUSE
white-light image, the S/N image and the segmentation map used to identify
continuum sources and regions of low S/N.  Figure 14: ORIGIN PCA thresholds.
The left image displays the field segmentation in 6 zones. The 6 plots show
the corresponding histograms of the PCA test statistics for each zone. The red
curves exhibit the Gaussian fit of the distribution. The computed thresholds
are displayed as blue vertical lines. They identify the minimum value above
which the spectrum triggers the PCA cleaning process. This value corresponds
to the 10% upper quantile of the empirical distribution.  Figure 15: Number of
iterations where each spaxel was involved in the iterative PCA is shown in the
left panel. The ORIGIN maxmap, computed as the maximum of the General Likehood
Ratio test datacube over the wavelength axis, is displayed in the right panel.
Figure 16: Counts of positive and negative local maxima in the GLR datacube as
a function of the threshold are shown in the left panel in blue and orange,
respectively. The GLR threshold value computed to achieve the target purity
value of 0.8 is displayed as a dotted vertical line. A similar plot for the
S/N datacube and a target purity of 0.90 is shown on the right panel.

###  5.3 ODHIN HST-prior spectra extraction

For deep exposures, source blending increases, to a point where source
confusion begins to be a problem. Source blending is already present in the
10-hour depth MOSAIC field and of course has an increased impact in the ten
times deeper MXDF data. Fortunately, the 3D information content of MUSE helps
to identify source mixing. For example, the location of the emission line can
be measured to an accuracy of about half a spaxel
($0\aas@@fstack{\prime\prime}1$) in reconstructed narrow bands when the S/N is
sufficiently high. For blended sources with different emission lines it is
then possible to separate the different contributions. Although this does not
work with the continuum, we have another very valuable piece of information:
the high spatial resolution broadband images from HST.

The combination of data at different spatial and spectral resolutions is known
as pansharpening in signal processing of hyper-spectral data (see the review
in Loncan et al., 2016). The aim of these methods is to reconstruct a datacube
with high spatial and spectral resolution. However, they are not directly
applicable to our problem due to the very different spatial and spectral
resolution in MUSE and HST101010There is a factor of $\approx 5$ in spatial
resolution and $\approx 1000$ in spectral resolution between MUSE and HST. and
the presence of noise. We have therefore developed a new and more robust
method, nicknamed ODHIN, based on linear source demixing.

The method assumes that at HST spatial resolution a source can be modeled as a
single spectrum. This is obviously not valid for nearby galaxies that exhibit
spatially resolved spectral variations. This is also not strictly valid at
high redshift where some line emission may be more extended than the broadband
morphology (e.g., Leclercq et al., 2017). ODHIN is similar to the TDOSE method
(Schmidt et al., 2019) developed for the MUSE-Wide survey (Urrutia et al.,
2019; Schmidt et al., 2021) but it differs in three aspects: it is
nonparametric, it uses multiple broadband HST images, and it implements a
regularization process to avoid noise amplification for very close sources.
The details of the algorithm are presented in Appendix C.

The inputs to the algorithm are the 3 HST broadband images (F606W, F775W, and
F850LP) that fall within the MUSE spectral range (Fig. 1), the HST source
catalog and its segmentation map and the MUSE datacube. We ran ODHIN for three
HST catalogs: our primary catalog R15 (Rafelski et al., 2015) but also the
CANDELS v2 (Whitaker et al., 2019) and 3D-HST (Skelton et al., 2014) catalogs.
Although these extractions use the same HST broadband images (XDF, Illingworth
et al., 2013), they are based on different segmentation maps and may thus
provide different results.

An example of source de-blending is given in Fig. 17. Two close sources can be
identified in the HST image (upper left panel) separated by 0.4 arcsec, but
they are blended at the MUSE spatial resolution as shown in the MUSE white
light image (second upper panel). The brightest source (F775W AB 27.2)
MID-7636 (RID-25333) is a nearby star-forming galaxy at $z=0.5$, with multiple
emission lines: for example, [O ii]3726,3729, [O iii]4959, H and [O
iii]5007(central row of Fig. 17). The second source (F775W AB 27.8) MID-6294
(RID-22260) is a $z=5.5$ $\mathrm{Ly}\alpha$ emitter with a strongly
asymmetric Ly line (last row of Fig. 17). The fact that the narrow-band image
of the Ly is slightly offset from the [O iii]5007 narrow-band image indicates
that the Ly emission is indeed related to RID-22260 and not RID-25333.
However, the overlap is too large to de-blend it from the narrowband
segmentation at the spatial resolution of MUSE. Therefore, for each source,
the spectra provided by the ORIGIN extraction (shown in orange in Fig. 17)
show strong contamination due to the nearby source. In contrast, the ODHIN
extraction (shown in blue) is able to de-blend the contributions from each
source.

One can see that, despite its intrinsic limitations111111We assume that the
broadband light distribution is representative of any light emission and do
not allow for spectral variation within a source., ODHIN is capable of de-
blending closed sources with good accuracy. We note, however, that in the
presence of emission that extends much beyond the HST broadband, the method
will miss part of the flux (see Fig. 19 for an example). Furthermore, the
method is by construction blind to any source undetected by HST.

Figure 17: Example of ODHIN de-blending of two close sources MID-7636
(RID-25333), a $z=0.5$ star forming galaxy, and MID-6294 (RID-22260), a
$z=5.5$ $\mathrm{Ly}\alpha$ emitter. The separation of the sources is 0.4
arcsec. The top row of panels shows, from left to right: the HST F775W, MUSE
white-light, [OIII]5007Å and Ly narrow band images ($\rm 5\times
5\,arcsec^{2}$ size). The central and bottom panels display the spectra
derived from ODHIN (in blue) and ORIGIN extraction (in orange) for each
source. Fluxes and observed wavelengths are expressed in units of
$\mathrm{10^{-20}\,erg\,s^{-1}\,cm^{-2}\,\AA^{-1}}$ and Å, respectively.

###  5.4 Redshift determination with pyMarZ

The software pyMarZ is a Python implementation of the redshift fitting
software MarZ (Hinton et al., 2016) originally developed for multi-object
fiber spectroscopy with OzDES (Lidman et al., 2020). MarZ itself is based on a
modified version of the cross-correlation algorithm AUTOZ (Baldry et al.,
2014). It is written in JavaScript with a dedicated user interface. The
development of the Python version was motivated by the need to integrate it
more easily into the full data analysis chain where all parts are written in
Python. We have also developed our own interface SourceInspector (Sect.
5.7.2), which is run separately from pyMarZ.

Compared to the version used in DR1 (Inami et al., 2017), the main change,
apart from the Python language, is the update of the spectral templates based
on a set of higher S/N spectra of selected sources in UDF-10 and MOSAIC
datacubes. We tested it on a selection of UDF-10 high confidence isolated
sources and found that the first ranked solution returned by pyMarZ is correct
84% of the time. We note that the probability of finding the correct solution
among the first two pyMarZ solutions is 95%. In the previous version used for
DR1, the probability of finding the correct solution was only 71% and 87% for
the first and first two returned solutions, respectively. Although the
improvement is significant, the selection of the correct redshift solution
requires additional information and cannot be performed automatically based
only on the top ranked pyMarZ solution. For this reason, we retain the top
five pyMarZ solutions for the manual inspections presented in Sect. 5.7. The
software, its documentation and the templates used for this work are made
publicly available (Appendix A.5).

###  5.5 Line flux measurements with pyPlatefit

In the past years, the increasing number of galaxy surveys with wide field or
highly multiplexed IFUs (e.g., SAMI, MANGA, and MUSE) has been accompanied by
the development of software tools to efficiently fit emission and/or
absorption lines in spectra or datacubes (e.g., pPXF Cappellari 2017; Sarzi et
al. 2017, ALFA Wesson 2016, LZIFU Ho et al. 2016, pipe3D Sánchez et al. 2021).
Most of these tools have been developed for the spatially resolved study of
nearby galaxies and are not optimized for the population of medium to high
redshift galaxies present in our observations. We have therefore developed a
specific fitting code pyPlatefit, inspired by the Platefit IDL code developed
for the SDSS survey (Tremonti et al., 2004; Brinchmann et al., 2004).

pyPlatefit is a python module to perform emission line fitting of astronomical
spectra. Although it was specifically developed for this data release121212We
used the Platefit IDL code for the DR1 release., it should work for spectra
delivered by other instruments, provided that the input spectrum is in MPDAF
format (Piqueras et al., 2019). The program takes as input a spectrum and an
input redshift and returns a list of emission line parameters: for example,
flux, velocity dispersion, equivalent width and their errors. The program
performs a continuum fit, using a simple stellar population model (Bruzual &
Charlot 2003, see Brinchmann et al. 2013), and then fits the emission lines
after subtracting the continuum. In addition, the program can also perform a
basic fit of the absorption lines. Lines are grouped into families, each of
which is fitted independently. All lines belonging to the same family share
the same velocity (or redshift) and velocity dispersion. The default lines
table (Appendix D) has 3 different families: the Balmer series, the non-Balmer
emission lines and the ISM absorption lines. In addition, some major resonant
lines like Ly or Mg ii are fitted separately as distinct families.

All emission and absorption lines are modeled as Gaussian with 3 parameters:
central wavelength ($\lambda_{0}$), wavelength dispersion ($\sigma$) and flux
($F_{0}$). The only exception is the Ly line, which is fitted as a skewed
Gaussian (Birnbaum, 1950), with an additional unit-less parameter ($\gamma$)
measuring the asymmetry of the line:

$F(\lambda)=F_{0}\left[1+\operatorname{erf}\left(\gamma\frac{\lambda-\lambda_{0}}{\sqrt{2}\sigma}\right)\right]\exp\left(-\frac{\left(\lambda-\lambda_{0}\right)^{2}}{2\sigma^{2}}\right).$
(4)

Optionally, a double peaked Ly line profile can be fitted as the sum of two
skewed Gaussians. The program uses the front-end lmfit minimization
package131313lmfit.github.io/lmfit-py (Newville et al., 2014) based on scipy
optimize nonlinear least-squares minimization (Virtanen et al., 2020). The
default minimization routine is the trust region reflective least squares
algorithm, which is more robust and gives a better error estimate than the
classical Levenberg–Marquardt algorithm. Other minimization algorithms can be
selected (see the lmfit documentation for details). The Markov chain Monte
Carlo emcee package (Foreman-Mackey et al., 2013) is also used to improve
error computation. Information about the public released version of the
software and its documentation is given in Sects. A.5 and D.

###  5.6 Narrowband measurements

For each redshift solution, we compute the narrowband images associated with
each detected emission or absorption line, as provided by pyPlatefit. These
narrowband images are important to confirm the detection of the line and to
verify its origin. For example, spectral pollution due to nearby sources will
be more easily spotted by looking at the location of the corresponding
narrowband image (e.g., Fig. 17). Another example is the case of a low S/N
spectrum when a noise pattern due to sky subtraction residuals is fitted by
pyPlatefit because its spectral location matches the expected wavelength of a
line. In this case, the narrowband image should not present a coherent peak at
the source location, thus helping discard the line measurement.

The narrowband process selects each emission and absorption line with
$\rm(S/N)_{line}>2$ plus all major141414e.g., Ly or [O ii]3726,3729; see
pyPlatefit input line table in Appendix 30. emission lines, regardless of
their S/N. The continuum is then subtracted from the datacube using a median
filter size of 101 Å total width (81 pixels). A narrowband image is finally
obtained using a simple summation of the continuum-subtracted datacube over
the wavelength range defined by the fitted center of the line and a spectral
width of $\rm\pm 2\times FWHM$, where FWHM is the fitted line width at half
maximum. For emission line doublets like [O ii]3726,3729, the wavelength range
is extended to cover both lines. For each narrowband image with
$\rm(S/N)_{line}>2$, we perform a SExtractor segmentation and derive its
center, flux and signal-to-noise ratio $\rm(S/N)_{nb}$.

We also derived composite emission and absorption narrowband images obtained
by optimally adding all emission or absorption narrow bands. The process
starts with the highest $\rm(S/N)_{nb}$ of the emission or absorption narrow
bands and adds additional narrowband images until the $\rm(S/N)_{nb}$ of the
composite narrowband stops increasing.

###  5.7 Visual inspections and reconciliations

The material produced in the previous steps, that is, for each candidate
source, the five redshift solutions provided by pyMarZ as well as their
emission and absorption lines fits and their associated narrow bands, must be
organized and evaluated by a group of experts. Once the conflicts have been
resolved during the reconciliation meetings, one can proceed to the creation
of the final catalog by assigning a unique identifier to the detected sources.
The details of the process are presented in the following sections. The
figures and statistics quoted are given for the MOSAIC largest data set, but
the process was similar for the UDF-10 and MXDF data sets.

####  5.7.1 Inspection sample

The MOSAIC parent sample is composed of 10450 sources divided into 7977 ODHIN
sources derived from the R15 catalog (Sect. 5.3) and 2473 ORIGIN sources
(Sect. 5.2). Previous investigation has shown that all DR1 ORIGIN sources to
which a secure redshift could be assigned have a purity greater than 0.9. We
therefore restrict the ORIGIN sample to the 1946 sources with purity greater
than 0.9.

The ORIGIN subsample contains all line sources that have at least one
detectable emission line, regardless of their continuum flux. The remaining
sources of interest are those whose continuum flux is bright enough to
identify other spectral features (e.g., absorption lines, breaks). After some
trial, we adopted a continuum S/N cut151515The continuum S/N is computed as
the maximum of the S/N spectrum after running a median filter of 100 pixels
width to remove the impact of emission lines. of $\mathrm{{0.8}\times
10^{-20}\,erg\,s^{-1}\,cm^{-2}\,\AA^{-1}}$ per spectral pixel, which
corresponds roughly to an AB F775W magnitude of 27.6. This reduces the ODHIN
sample to 782 sources161616Our experience has taught us that galaxies with
fainter continuum cannot be attributed an aborption line redshift and that
galaxies with emission lines are already detected with ORIGIN..

We further reduced the MOSAIC inspection sample by removing all sources
located in the UDF-10 area, taking a margin of 1 arcsec. The final MOSAIC
inspection sample is composed of 2412 sources: 1727 ORIGIN and 685 ODHIN
sources. This sample was randomly distributed to four groups of three experts.
In addition, a randomly selected common subsample of 50 sources was given to
the four groups. This control sample will be used to assess the homogeneity of
the classification between the groups. We proceeded likewise for the UDF-10
and MXDF data sets and produced final (i.e., ODHIN, ORIGIN, and control)
inspection samples of 1234 and 1599 sources, respectively.

####  5.7.2 SourceInspector

The inspection package provided to each expert consists of a small database, a
set of interactive html files and the SourceInspector software.
SourceInspector is a PyQt tool that allows the expert to display each
candidate source, to select the redshift solution or to provide a new one, to
give a confidence to the redshift assignment and to match the source to other
sources (e.g., ORIGIN to ODHIN) and to HST catalogs. All operations are
performed locally. An example of the SourceInspector interface is shown in
Fig. 18. The inspection results are exported at the end of the process to be
used later for the reconciliation step.

Figure 18: SourceInspector interface in evaluation mode. The left side of the
interface displays the interactive html file for the current source. It shows
the five redshift solutions proposed by pyMarZ plus their corresponding
pyPlatefit line fits and associated narrow bands in the bottom panels. The
html upper panels display MUSE and HST images, and the locations of
neighboring sources from HST catalogs. The right side of the interface allows
the user to navigate in the source list, select the redshift solution or
assign a new redshift, assign confidence and provide crossmatch information.

####  5.7.3 Redshift confidence

The expert must assign a confidence to the redshift solution. This redshift
confidence (called ZCONF in the tables) is defined as follows.

A confidence level of 0 means no redshift solution can be found. It should be
noted that this is decoupled from the detection status. An emission line can
be well detected with high purity by ORIGIN, even though no redshift solution
can be found.

A confidence level of 1 indicates a low confidence. A tentative redshift
solution has been found. Reasons for assigning $\rm ZCONF=1$ include low S/N
of the lines (e.g., $\rm S/N_{line}<3$), poor fit, existence of other valid
redshift solutions, noisy associated narrow bands, additional lines present in
the spectrum with reasonable S/N but unexplained by the proposed redshift
solution.

A confidence level of 2 indicates a good confidence. For non
$\mathrm{Ly}\alpha$ emitters (i.e., $z<2.9$) sources, we assign this
confidence level to sources with multiple lines detected with good S/N ($\rm
S/N_{line}>5$). For example, a resolved [O ii]3726,3729 doublet with a clear
narrowband would be sufficient to assign $\rm ZCONF=2$. For
$\mathrm{Ly}\alpha$ emitters, the requirement is to have a Ly line with good
S/N (i.e., $\rm S/N_{line}>5$) and a width and asymmetry compatible with Ly
line shapes.

Finally, a confidence level of 3 indicates a high confidence. For
non-$\mathrm{Ly}\alpha$ emitters, the criteria are similar to the $\rm
ZCONF=2$ requirements, but we expect more lines and higher S/N, as well as
high S/N narrowband images. For $\mathrm{Ly}\alpha$ emitters, if there is no
other line than Ly, we require the Ly line to achieve high S/N ($\rm
S/N_{line}>7$) with the expected line shape: that is a pronounced red
asymmetrical line profile171717$\gamma$ is the asymmetric Gaussian parameter
used by pyPlatefit (Sect 5.5). ($\rm\gamma>2$) and/or a blue bump or double
peaked line profile181818Double peaked Ly profiles are easily differentiated
from [O ii]3726,3729 or C iii]1907,1909 using the line spacings and widths..

In the case of ORIGIN detections and if the source can be matched to an HST
counterpart, it adds confidence to the detection. If, in addition, the R15
photometric redshift is reliable and matches well with the MUSE redshift (see
Eq. 5), this adds confidence to the redshift assignment and could lead to a
higher ZCONF than in the case where the ORIGIN source is not matched with an
HST source. On the other hand, photometric redshifts for faint sources are
often not very reliable and thus disagreement between the MUSE and photometric
redshifts is not considered a strong negative constraint.

We note that in the case of low S/N spectra with a single emission line, we
are not always able to determine the nature of that line with a high degree of
certainty. In many cases, the combination of line shape and photometric
redshift argues in favor of Ly, but it is desirable to have a handle on the
amount of misclassifications in this subsample. We focus on [O ii] and Ly
emitters here as other emission lines are much less likely because of the
small accessible volume.

In order to estimate the expected fraction of Ly and [O ii] emitters as a
function of magnitude, we need LFs and rest equivalent-width distribution
functions (EWDFs) for both. We take the Ly LF from Herenz et al. (2019) and
the Ly EWDF from Hashimoto et al. (2017). For [O ii] we adopted the evolving
LF from Saito et al. (2020), which agrees fairly well with our UDF data. We
also need an EWDF for [O ii], which we prefer to derive from our high quality
($\rm ZCONF\geq 2$) data. To do that we adopt the point-source selection
function from Herenz et al. (2019), which we derived for Ly emitters but which
reproduces the minimum [O ii] flux as a function of wavelength fairly well for
the present sample when the selection function is adjusted for difference in
exposure time by $\rm\sqrt{t}$. We then use a simple Vmax estimator for the
EWDF after verifying that the [O ii] line luminosity and equivalent width are
approximately independent, in good agreement with previous work (Blanton &
Lin, 2000; Ciardullo et al., 2013). The distribution is found to be well
fitted by a Gaussian in log EW(Å) with mean = 1.61 and sigma = 0.17.

Taking these together with the assumption that the spectrum of our galaxies is
approximately flat in flux density ($\rm f_{\nu}$) we can calculate the
expected number of [O ii] and Ly emitters with observed magnitude F775W
fainter than a particular value. When we do that, we find that at magnitudes
$\rm F775W>28.5$ the number of $\mathrm{Ly}\alpha$ emitters is expected to be
$10$—$100$ times that of the [O ii] emitters at all redshifts of relevance
down to a flux limit of $\mathrm{10^{-18}\,erg\,s^{-1}\,cm^{-2}\,\AA^{-1}}$.
Thus, we conclude that the likelihood of misclassified lines is less than 1/10
for galaxies fainter than 28.5.

####  5.7.4 Cross identifications

For ORIGIN sources, the expert must also identify the possible HST
counterparts and, likewise, must find the possible ORIGIN counterparts for
ODHIN sources. Crossmatching ORIGIN and HST sources can be difficult when the
field is overcrowded or if the ORIGIN location is offset with respect to HST
(see an example in Fig. 37). The use of the photometric redshift can help with
the HST association, but not always. At a later stage we perform an evaluation
of the HST crossmatch confidence (Sect. 5.8.2).

Occasionally a source is visually identified in one of the HST images at the
ORIGIN location, but is not present in the HST R15 catalog. This might occur
when the source has been incorrectly merged with a nearby source by
SExtractor, but in other cases the source is simply too faint and escapes
detection. We identify these cases by setting ”missed, faint, or undetect” in
the corresponding HST crossmatch flag (see Sect. 6.2).

The ORIGIN process identifies individual emission lines and merges them into a
single source when they are close (see Sect. 3.5 in Mary et al. 2020).
However, sometimes the process fails to merge emission lines when the source
is very extended, or, on the contrary, it merges lines that belong to two
different sources (i.e., superpositions). The expert is then responsible for
evaluating if the ORIGIN sources should be split or merged with another ORIGIN
source. Such actions are then performed in the final step (Sect.5.8).

Cross-identification with the DR1 catalog is also performed during this
evaluation process. When a source is found at the location of a DR1 source but
is assigned a different redshift, we mark it as a DR1 unmatched source. On the
other hand, if only the redshift confidence is different, we mark the source
as DR1 matched.

####  5.7.5 Reconciliation

When all experts in the group have completed their evaluation, the group
leader performs the reconciliation process. All sources that disagree on any
of the assigned items (i.e., redshift value, redshift confidence, HST, DR1, or
ORIGIN crossmatch) are flagged for the reconciliation meeting. In addition,
any source flagged with ”revisit” by an expert is added to the review list.
The only exception to this process is when the only discrepancy is due to
redshift confidences between two and three, in which case the source is not
flagged for reconciliation and its confidence is set to the average confidence
rounded to the nearest integer.

An example of the statistics returned by SourceInspector in automatic checking
mode is given in Table. 5. Disagreements are then discussed and resolved in a
face-to-face reconciliation meeting with the help of the reconciliation mode
of SourceInspector. On average, between 30% and 40% of the original sample had
to be revisited.

Most of the sources to be re-examined are weak, single-line sources for which
experts have a different assessment of the reality of the measured redshift.
For example, some less experienced experts would be inclined to identify a
$\mathrm{Ly}\alpha$ emitter for a low S/N line, even if there is no additional
information to support this assumption. In most cases, however, when
reconciling, the redshift solution will be ruled out.

Table 5: Example of automatic reconciliation.

MOSAIC group 1 reconciliation |
---|---
Total number of sources | 790
Sources autochecked with $\rm ZCONF=0$ | 57
Sources autochecked with $\rm ZCONF>0$ | 430
Sources to be manually reconcilied | 303
With revisit flag | 24
With discrepant redshift | 185
With discrepant redshift confidence | 45
With discrepant matching information | 36
With discrepant split or merge info | 13
191919 ”autochecked” means that all experts agree.

####  5.7.6 Control sample

A control sample of 50 randomly selected MOSAIC sources was given to the same
four experts groups. The same experiment was done for the MXDF evaluation.
After reconciliation, the control sample was analyzed for remaining
disagreements. The 50 sources consisted of 36 ORIGIN and 14 ODHIN sources. All
groups agreed on the redshift assessments. Eight sources did not obtain a
redshift solution. Of the 42 sources with a redshift solution, nine sources
were credited with different redshift confidence: five sources were ranked
with low and good/high confidence and four sources were ranked with good and
high confidences. Four sources were also crossmatched to different HST
sources.

The results for this control sample demonstrate that we achieved overall
excellent agreement among all experts, with no discrepancies in redshift
solutions. Redshift confidence is more subjective and we do identify 20% of
disagreements between groups. We note that we do not consider the difference
between confidence two and three to be critical. Nevertheless, we have 10% of
sources that were considered low confidence by some groups and good or high
confidence by others.

Similar results were obtained on a control sample of 49 MXDF sources given to
three different groups. All groups agreed on the redshifts, but there were ten
sources with different confidence assignments. In addition, some groups
assigned a low confidence level to 10 sources, while the others classified
them as having no redshift (confidence 0).

As we discuss in the next section, a review of all sources with redshifts is
performed as a final step by a single expert. This helps to homogenize the
results, but given the amount and complexity of information that goes into the
confidence evaluation (Sect. 5.7.3), it currently seems difficult to have a
more objective criterion for redshift confidence assignment.

###  5.8 Creation of final data products

In this final step, we create the final catalog by giving a unique MUSE
identifier to all sources with an assigned redshift. We start by splitting the
catalog in ORIGIN only, ODHIN only and crossmatched ORIGIN-ODHIN sources. All
sources already matched to the DR1 catalog keep the same MUSE identifier, all
others are given a unique new identifier.

Split and merge of ORIGIN sources defined during the inspection step was
implemented. A final review of all sources with redshift was then completed.
During this review the reference extraction (Sect. 5.8.1) and the reference
center (Sect. 5.6) were chosen. Emission or absorption line fits were fine-
tuned with pyPlatefit when needed and double asymmetric Ly fits implemented
for double peaked profiles. The computed narrowband $\rm S/N_{nb}$ was used as
an additional guide for the redshift confidence (Fig. 23), which was updated
when requested.

####  5.8.1 Selection of reference extraction

For each source, we derive the spectra from the MUSE datacube using different
extraction schemes (e.g., ORIGIN or ODHIN). There is no single way to extract
the signal: for example, the stellar continuum is not necessarily superimposed
on the location of the flux peak of the emission lines, and its spatial extend
is often different. The most obvious cases are those of $\mathrm{Ly}\alpha$
emitters, which show extended Ly emission that is on average ten times larger
than their continuum counterparts (e.g., Wisotzki et al. 2016; Leclercq et al.
2017). In this case, a spectrum based on the continuum surface brightness will
miss a large faction of the Ly flux, and vice versa, a spectrum optimally
extracted for the Ly emission will be suboptimal for the continuum S/N.

In the ORIGIN extraction scheme (Sect. 5.2), the pseudo-narrowband image is
derived from an improved continuum subtracted datacube using an iterative
principal component (PCA) scheme and thus has the advantage of being mostly
free from continuum residuals. It also has a high S/N, as it results from a
process similar to optimal filtering (see Sect. 3.2 of Mary et al. 2020).
However, in the case of a very bright emission line, the PCA process may
produce artifacts by removing part of emission line. In this case, we perform
an optimal extraction based on the narrowband composite image and its
segmentation (Sect. 5.6). This narrowband extraction method is called NBEXT in
the rest of the document.

We illustrate the difference between the two methods with two examples in Fig.
19. Figure 19(a) shows the results of the ORIGIN and NBEXT weighted
extractions for a bright $\mathrm{Ly}\alpha$ emitter in which this effect is
particularly strong. It can be seen that in the case of the ORIGIN pseudo
narrowband image (upper central panel), the bright central Ly emission has
been partially removed by the PCA, leaving only the outskirts. As a result,
only a small portion of the Ly flux is recovered by the weighted extraction
(shown in red in the lower panels). We note that the ODHIN extraction (Sect.
5.3) based on the R15 segmentation map (blue line in bottom panels) is more
efficient but still misses some Ly flux compared to the narrow band Ly
extraction (green line in bottom panels). The second example shown in Fig.
19(b) presents the same extraction methods applied to a faint
$\mathrm{Ly}\alpha$ emitter. We observe that the narrowband image is too noisy
to capture the Ly flux and that the ORIGIN weighted extraction performs
better. We stress that in most cases, however, the different extractions
methods give very similar results.

For each source, we assigned a reference extraction. In general, we favored
ORIGIN or NBEXT for $\mathrm{Ly}\alpha$ emitters and ODHIN for low-$z$ sources
when the continuum was not too faint. In some cases where the contamination is
strong, ODHIN was preferred for $\mathrm{Ly}\alpha$ emitters, even if it
misses some Ly flux. The selection was made on a case-by-case basis. However,
the reference extraction may not be the best for a given science case and so
we provide all extractions in the delivered data products (Appendix A.3).

(a) MID-106.

(b) MID-7295.

Figure 19: Example of sources where spectral extraction gives different
results. Three different extraction schemes are shown for the
$\mathrm{Ly}\alpha$ emitters MID-106 (a) and MID-7295 (b). For each subfigure
we show the full spectrum (lower left panel) and a zoomed region around the Ly
line (lower right panel). In green we display the spectrum derived from the
weighted extraction using the narrowband Ly weighted map (NBEXT) and its
segmentation (top left panel), in red the ORIGIN weighted extraction (top
central panel) and in blue the ODHIN extraction based on the R15 segmentation
map (upper right panel). The colored area in the left and central top panels
indicate the source (orange) and masked objects (magenta) segmentation maps.
In the case of the low S/N source MID-7295 (b), the spectra have been smoothed
with a 5 pixels width kernel and a Gaussian filter of 0.4 arcsec FWHM has been
applied to the Ly narrowband image (top left panel of subfigure b). Observed
wavelengths are in Å units.

####  5.8.2 Crossmatch with HST source catalogs

During the inspection step (Sect. 5.7) the ORIGIN sources were matched to the
HST R15 catalog. However, since the publication of Rafelski et al. (2015), new
HST source catalogs in the CANDELS area have been published. These catalogs
are based on an improved data reduction of HST photometry (Whitaker et al.,
2019), additional HST GRISM information as in 3D-HST (Skelton et al., 2014),
or additional ground-based photometric bands (Merlin et al., 2021). In
addition, the detection image202020The Rafelski et al. (2015) detection image
is the average of the four optical ACS images with the four WFC3 images.
Whitaker et al. (2019) perform a noise-equalized combination of one ACS
(F850LP) with three WFC3 bands (F125W, F140W, F160W). The Merlin et al. (2021)
catalog is based on the Guo et al. (2013) catalog, which uses the F160W as
detection image. and the segmentation input parameters were different,
resulting in different segmentation maps. The biggest difference is for
closely packed faint galaxies, which can be merged in a single source in some
catalogs, whereas they are split into subcomponents in others. In our
investigation, we found many cases where the HST source should be split into
two different sources (e.g., Fig. 28) or, on the contrary, where two close HST
sources should be considered as a single one.

In this final step, we therefore crossmatched our results with the CANDELS v2
(Whitaker et al., 2019), 3D-HST (Skelton et al., 2014) and ASTRODEEP (Merlin
et al., 2021) catalogs. These catalogs and their corresponding segmentation
maps were also used to produce the corresponding ODHIN extraction (Sect. 5.3).

For all sources, we calculated an HST environment index to quantify whether
the source can be considered spatially isolated or has close or distant
neighbors. This index is the percentage of the area covered by neighboring
sources measured in the R15 segmentation map212121We use the CANDELS v2
segmentation map when the HST source is not present in R15. in a given
aperture centered on the MUSE source location. We used aperture radii of 0.5
and 1.0 arcsec to define the close and distant environment index,
respectively. At the MUSE spatial resolution, an object with a nonzero close
environment index will overlap with its neighbors, at least partially, while
an object with zero close environment index and nonzero distant environment
index will be resolved by MUSE. A source is considered isolated if its distant
environment index is zero.

The presence or absence of close or distant neighbors is coded in the IFLAG
column of the main catalog (Sect. A.2). It is defined as follows: (i) Sources
with a close HST environment index of zero and a distant environment index of
less than 5% are considered isolated with respect to HST. If, in addition,
there is no other non-HST DR2 sources within 1 arcsec radius, they are defined
as fully isolated. Their isolation indicator ICONF is set to three. (ii)
Sources with $\rm ICONF=2$ have only distant neighbors, either from HST or
from DR2 non-HST sources. (iii) All sources with at least one close HST or
non-HST neighbor have ICONF set to one.

To assess the quality of the provided matching information, we assigned a
matching confidence indicator (hereafter MCONF) to each source. Its values are
given in Table 6 and range from 0 (ambiguous) to 3 (high confidence). MCONF
combines the neighborhood information as defined above with additional
redshift and narrowband offset information.

An HST source is matched in redshift with a MUSE source when its photometric
redshift is considered reliable and the MUSE redshift is within the 95%
confidence interval, that is,

$\frac{z_{max}-z_{min}}{1+z_{best}}<0.22\\\ z_{min}<z_{MUSE}<z_{max},$ (5)

where $z_{best},z_{min},z_{max}$ are respectively the photometric redshift
solution and its bounds given by the 95% error probability. In the case of
3D-HST spectroscopic redshifts we use $|z_{MUSE}-z_{spec}|/(1+z_{MUSE})<0.01$.

The last criterion is the offset between the MUSE absorption or emission
composite narrowband (Sect. 5.6) from the center of the matched HST source.
This offset is compared to the size of the HST source and is considered
acceptable when it satisfies the following condition: $\rm
offset<\max{(2\times size,0.4\arcsec)}$ where size is the radius of the
circular HST aperture (in arcsec) enclosing half of the total flux, as given
in the catalog.

All isolated sources that do not have an HST counterpart in any catalog are
assigned a high matching confidence ($\rm MCONF=3$). There are 99 sources in
this category. All isolated sources matched to a single HST counterpart are
also assigned a MCONF value of three if all reliable photometric redshifts of
the two used HST catalogs can be matched to the MUSE redshift and if there is
no significant offset between the HST source and the narrowband location. If
the offset is too large and/or there is at least one reliable photo-z
discrepant with respect to MUSE, MCONF is set to two.

When sources cannot be considered isolated, in addition to the previous
criteria, we search if any neighbor, not already assigned to any MUSE sources,
has a reliable photo-z that could be matched to the MUSE redshift. If so,
MCONF is lowered to two or one, depending on if the neighboring source is
distant or close to the MUSE source. Finally, sources with multiple HST
counterparts assigned during the manual evaluation are assigned a zero value
of MCONF. There are only eight sources in this category.

The counts of the matching confidence are given in Table. 6. Most MUSE sources
have a high HST matching confidence (68%), either because the MUSE source can
be assigned to a unique identifier in the HST catalogs (53%) or because they
are, instead, isolated and have no HST counterpart in the catalog (15%).

Table 6: Isolation flag (IFLAG) and HST crossmatch confidence (MCONF) definition and counts. IFLAG | Description | Counts | MCONF | Description | Count
---|---|---|---|---|---
|  |  | 0 | Ambiguous | 8
1 | Close neighbors | 679 | 1 | Low confidence | 61
2 | Distant neighbors | 154 | 2 | Good confidence | 644
3 | Isolated | 1388 | 3 | High confidence | 1508

####  5.8.3 Systemic redshift estimation for $\mathrm{Ly}\alpha$ emitters

For most $\mathrm{Ly}\alpha$ emitters the measured redshift is based on the
peak of the Ly line. However, due to the resonant scattering properties of Ly
photons in the interstellar medium, the Ly redshift is systematically
different from the systemic redshift (e.g., Shapley et al., 2003; McLinden et
al., 2011; Rakic et al., 2011; Song et al., 2014). Typical velocity offsets
are $\approx$ 200 $\mathrm{km\,s^{-1}}$ for $\mathrm{Ly}\alpha$ emitters with
larger values ($\approx$ 500 $\mathrm{km\,s^{-1}}$) for Lyman-break galaxies
(Shibuya et al., 2014; Muzahid et al., 2020).

We used the empirical recipes provided by Verhamme et al. (2018) to estimate
the velocity offset between the Ly and systemic redshift. In their paper,
Verhamme et al. (2018) supply two types of corrections. The first correction
is based on the separation between the two peaks of the Ly line (i.e., the
blue bump separation). The center of the two peaks is a good estimate of the
systemic redshift with a rest-frame scatter of 53 $\mathrm{km\,s^{-1}}$. The
asymmetric double peaked line profile used in pyPlatefit (Sect. 5.5) provides
this value directly. We simply add the 53 $\mathrm{km\,s^{-1}}$ standard
deviation derived by Verhamme et al. (2018) to the redshift error budget in
quadrature.

The second type of correction is based on the FWHM of the Ly line profile. The
rest-frame velocity offset ($\rm\Delta V$ in $\mathrm{km\,s^{-1}}$) is
inferred from the following linear relation with the rest frame measured FWHM
in $\mathrm{km\,s^{-1}}$ (uncorrected for the MUSE LSF): $\rm\Delta
V=0.9\times FWHM-34\,\mathrm{km\,s^{-1}}$. It has a larger scatter of 72
$\mathrm{km\,s^{-1}}$ and therefore we use the first correction when a double
Ly peak was detected and successfully fitted.

We checked the accuracy of the correction by selecting the few Ly emitters
with high redshift confidence, high S/N Ly line and additional high S/N
nonresonant emission or absorption lines. This selection resulted in a sample
of 14 $\mathrm{Ly}\alpha$ emitters: six with single Ly line profiles and eight
with double peaked profiles. Systemic velocities were measured from forbidden
emission lines (mostly C iii]1907,1909) or strong absorption lines (e.g.,
Boogaard et al. 2021). The comparison of the velocity offsets before and after
correction is presented in Fig. 20. We can see that the systematic velocity
offsets of about 300 $\mathrm{km\,s^{-1}}$ are reduced to a few tens of
$\mathrm{km\,s^{-1}}$. The measured scatter (77 $\mathrm{km\,s^{-1}}$ RMS) is
comparable to the values measured by Verhamme et al. (2018).

Figure 20: Check of systemic redshift correction for a subsample of DR2
$\mathrm{Ly}\alpha$ emitters with systemic redshifts measurements. The
histogram of the rest-frame velocity difference between the measured Ly red
peak location ($\rm V_{Ly\textalpha}$) and the systemic velocity ($\rm
V_{sys}$) is shown in orange. The same histogram after systemic velocity
correction ($\rm V_{Ly\textalpha}^{*}$) is shown in blue. The type of
correction, that is, blue bump separation (bb) or Ly line FWHM (fwhm), are
identified in respectively darker and lighter colors. The dashed vertical line
displays the median value of the velocity offset distribution after correction
($\rm-47\,\mathrm{km\,s^{-1}}$).

####  5.8.4 Final catalogs

In Table 7 we display the source count by extraction method and data set. We
find that the majority of sources (67%, 1499/2221) detected with ORIGIN have
also an HST counterpart, but we also observe that the fraction of ORIGIN
sources without HST counterparts increases with depth, from 13% (202/1505) in
the MOSAIC to 29% (203/691) in the MXDF data sets.

When a source is detected in multiple data sets (e.g., UDF-10 and MXDF), we
retain the data set with the greatest depth at the source location. The only
exception is MID-171, where we selected the UDF-10 instead of MXDF data set,
despite its lower depth (31 versus 136 hours). The reason for this is the Ly
line that falls into the MXDF masked wavelength range because of the sodium
laser guide star notch filter (Fig. 1). When available, the other data sets
are also recorded in the delivered data product (Appendix A.3).

The final catalog contains 2221 sources with redshifts. The statistics and
properties of these sources are discussed in the following sections.

Table 7: Source counts with respect to detection method.

data set | ORIaa$a$ORIGIN | ODHbb$b$ODHIN |  ORI & ODH | Total
---|---|---|---|---
MOSAIC | 202 | 220 | 1083 | 1505
UDF-10 | 78 | 48 | 225 | 351
MXDF | 203 | 78 | 410 | 691
DR2 | 428 | 294 | 1499 | 2221
222222

##  6 Results

###  6.1 Redshifts

One of the most visible results of the survey is the number of galaxies for
which we were able to assign high quality redshifts. A detailed presentation
is given in the following subsections.

####  6.1.1 Redshift distribution and confidence

The distribution of redshifts is shown in Fig. 21. Table 8 presents the
breakdown of the 2221 redshifts in five categories. Given the small area (8.3
arcmin2) and the Galactic pole location, we identified only eight stars,
mostly M types. Likewise, the number of nearby galaxies ($\rm z<0.25$) is
limited (25) in this small volume (280 $\mathrm{cMpc^{3}}$). For this latter
category, dedicated wide field spectroscopic surveys (e.g., SDSS Blanton et
al. 2017 and MANGA Bundy et al. 2015) are much more efficient. Nevertheless,
for these few bright galaxies (median AB F775W of 23.3), we obtain spatially
resolved spectroscopy with exquisite S/N up to large galactocentric distance.
An example is given in Fig. 22 row a.

In the redshift range 0.25–1.5, the [O ii]3726,3729 doublet falls in the MUSE
wavelength range. At the MUSE spectral resolution ($\rm R\approx 3000$), the
[O ii] doublet is resolved in most cases, making line identification easier.
With 677 galaxies, this population is the second most populated category,
after the $\mathrm{Ly}\alpha$ emitters. These galaxies cover a broad range in
F775W AB magnitude with 23.3, 25.7 and 27.6 for the 10%, 50% and 90%
percentiles, respectively. Examples are given in Fig. 22, rows b and c. We
note that a small fraction (5%) of the galaxies in this wavelength range are
not star forming [O ii] emitters, but passive galaxies detected through their
absorption lines.

The redshift range 1.5–2.8 corresponds to the MUSE redshift desert, a region
where [O ii] is redshifted out of the MUSE wavelength range and Ly is not yet
observable. The number of galaxies drops to 201. In this redshift range,
galaxies are detected either through their UV absorption lines (50%) or by
fainter emission lines, mainly the C iii]1907,1909 or Mg ii 2797,2803 doublets
(see examples in Fig. 22, rows d and e).

The last category comprises the $\mathrm{Ly}\alpha$ emitters, covering the
redshift range 2.8–6.7. With 1308 $\mathrm{Ly}\alpha$ emitters or 60% of the
complete sample, it is the most numerous and also the faintest population
(29.4 median AB F775W). The vast majority of these galaxies (98%) are
identified by their high equivalent width asymmetric Ly emission line but a
few galaxies (3%) exhibit in addition fainter emission lines such as C iv
1548,1550 or C iii]1907,1909. A small percentage (2%) are Ly absorbers with
broad absorption and little or no detectable Ly emission. We note that a
significant fraction of the identified $\mathrm{Ly}\alpha$ emitters have no
HST counterpart (see Sect 6.2). A few examples of $\mathrm{Ly}\alpha$ emitters
are displayed in Fig. 22, rows f to h.

It can be seen in Table 8 and in the right panel of Fig. 21 that the fraction
of good to high redshift confidence ($\rm ZCONF>1$) is high (87%) for the low
and intermediate redshift populations (i.e., $z<2.8$). This ratio is slightly
lower (70%) for the $\mathrm{Ly}\alpha$ emitters for the reason that in most
cases the Ly line is the only line identified.

In Sect. 5.7.3 we presented the redshift confidence definition and assessment.
Although it is a complex process based on multiple parameters and data sets,
including the expert’s judgement, emission line S/N is one of the primary
factors in the decision process. This can be seen in Fig. 23 where we show the
confidence values as a function of the main emission line S/N derived from the
extracted spectra and the corresponding narrowband S/N.

We point out that sources with low redshift confidence ($\rm ZCONF=1$) have
nevertheless a clear detection signal, but in general the S/N is too low to
assign a redshift with high confidence. An example is given in Fig. 24 where
the same source is shown for the UDF-10 and MXDF data sets with 34- and
128-hour depths, respectively. In the UDF-10 data set (upper left panel in
Fig. 24), the emission line is detected with an S/N of 5.3 and does not appear
to be very asymmetric ($\rm\gamma=1.5\pm 2.6$). The corresponding narrowband
(lower left panel) shows a clear peak, although with a low S/N of 3.5. The
measured redshift is 3.17, which is not far from the R15 photometric redshift
of 2.74 estimated for the matched HST source. We note, however, that the
photometric redshift ($\rm zp$) measurement is not very reliable with
$\rm(zp_{max}-zp_{min})/(1+zp)=0.27$. Taking all these points into
consideration, the expert gave a low redshift confidence to the source. The
same source with the deepest MXDF observations (right column in Fig. 24)
obtains a higher S/N of 10.7 and 7.8 for the emission line and the
corresponding narrowband, respectively. The line profile is now clearly
asymmetric with $\rm\gamma=3.6\pm 2.6$. Given these measurements, a high
redshift confidence ($\rm ZCONF=3$) was given to the source.

Table 8: Distribution of redshifts by confidence and for the different classes
of sources.

Class | ZCONF
---|---
| 1 | 2 | 3 | 1-3
Stars | 0 | 0 | 8 | 8 | 0.4%
Nearby ($\rm z<0.25$) | 4 | 0 | 21 | 25 | 1.1%
[O ii] ($\rm 0.25<z<1.5$) | 55 | 122 | 500 | 677 | 30.5%
Desert ($\rm 1.5<z<2.8$) | 55 | 57 | 89 | 201 | 9.0%
Ly ($\rm z>2.8$) | 396 | 538 | 374 | 1308 | 58.9%
All | 510 | 717 | 994 | 2221 | 100.0%
232323 ZCONF range from 1 (less secure redshift) to 3 (most secure redshift);
see Sect. 5.7.3.

Figure 21: Redshift distribution of the MUSE data release II. The magnitude -
redshift scatter plot is displayed on the left panel. All MUSE sources with an
HST counterpart in at least one of the HST catalogs are shown as blue symbols.
We note that objects with low S/N in F775W could still be in HST catalogs if
they are detected in another filter (e.g., NIR filters for high-z objects).
The matched spectroscopic redshifts previously published (Sect. 6.1.3) are
shown as yellow circles. For sources without catalog HST counterparts, the
brown circles show the 0$\aas@@fstack{\prime\prime}$4 aperture magnitude as
measured on the HST F775 image (c.f Sect. 6.2). The horizontal dashed line
indicates the median F775W $\rm 3\sigma$ limiting magnitude. The redshift
histogram, colored by redshift confidence, is presented in the right panel.
Figure 22: Examples of sources. From left to right the columns display (1) the
HST image in the F775W or F850LP filter (2) the emission or absorption
narrowband for the main line (3) the spectrum zoomed-in on the main line (in
blue) and its fit (in red) (4) the full spectrum (in blue). The green curves
show the noise standard deviation. Except for the first two sources, the
narrowband images and the full spectra have been filtered with respectively a
Gaussian of 0.5″ FWHM and a box filter of 5 pixels. Image sizes are
$5\arcsec\times 5\arcsec$. Observed wavelengths are in Å and the flux unit is
$\mathrm{10^{-20}\,erg\,s^{-1}\,cm^{-2}\,\AA^{-1}}$. For each source the MUSE
identifier, depth in hours, redshift value and confidence, and AB HST
magnitude are indicated.  Figure 23: S/N estimates as a function of redshift
confidence. The log of the S/N emission or absorption composite narrowband is
shown as a function of the log of the leading line S/N derived from
pyPlatefit. Low (1), good (2) and high (3) redshift confidences are shown in
respectively orange, blue and green colors, while low-z ($z<2.9$) and high-z
($z>2.9$) sources are displayed in the left and right panels, respectively.
Figure 24: Impact of depth on the redshift confidence assignment. The Ly line
(data in blue, fit in red) and the corresponding narrowband image of the
$\mathrm{Ly}\alpha$ emitter MID-690 are shown for the UDF-10 (left column) and
the MXDF (right column) data sets. The green curves display the noise standard
deviation. The narrow bands have been smoothed with a Gaussian of 0.5″ FWHM.
Image sizes are $5\arcsec\times 5\arcsec$. Observed wavelength and flux units
are respectively Å and $\mathrm{10^{-20}\,erg\,s^{-1}\,cm^{-2}\,\AA^{-1}}$.

####  6.1.2 Comparison with DR1

We recall that the DR1 catalog (Inami et al., 2017) was based on an early data
reduction of the UDF-10 and MOSAIC data sets and did not include the new MXDF
data set. Furthermore, although the overall detection and inspection process
was similar, there are a few important differences. First, the data reduction
was significantly improved with the use of super-flatfield, improved sky-
subtraction and self-calibration (Appendix B). Second, the HST basic aperture
extraction was replaced by the de-blending ODHIN procedure (Sect. 5.3). Third,
we use new templates for pyMarZ (Sect. 5.4) and the improved line fitting
software pyPlatefit (Sect. 5.5). So we expect some differences between DR1 and
the new catalog, even in the area outside the MXDF footprint. A summary of the
differences is given in Tables 9(b) and 9(b).

The DR1 catalog contains 1574 sources. While most of the DR1 sources are
present in the new catalog, a small fraction (9%) were not confirmed. The
reasons for the failure are given in Tab. 9(b). The rejected sample is mainly
composed of low-confidence DR1 sources that were rejected or assigned a
different redshift in the new deeper (MXDF) or improved (MOSAIC, UDF-10) data
sets. For traceability purposes, new MUSE identifiers were given to DR1
sources that had been assigned a different redshift.

In Fig. 25 we display the quantitative improvement between DR1 and DR2. As
shown in Tab. 9(b), the majority of the 792 new DR2 sources (53%) are provided
by the MXDF data set. This was expected given the improved depth, that is, a
factor of 12 and 4 with respect to MOSAIC and UDF-10, respectively. However,
there is still a significant gain of about 25% for the UDF-10 and MOSAIC data
sets. This demonstrates that the effort spent on improving the entire process
(i.e., data reduction, source detection and classification) has paid off.

Table 9: Comparison of DR2 and DR1 catalogs. (a) DR1 sources not matched in
DR2

Reason | $\rm ZCONF=1$ | $\rm ZCONF>1$ | Total
---|---|---|---
New Z | 25 | 1.6% | 5 | 0.3% | 30 | 1.9%
No Z | 88 | 5.6% | 17 | 1.1% | 105 | 6.7%
Duplicate | 1 | 0.1% | 5 | 0.3% | 6 | 0.4%
New Match | 1 | 0.1% | 0 | 0.0% | 1 | 0.1%
All | 115 | 7.3% | 27 | 1.7% | 142 | 9.0%

(b) New DR2 sources with respect to DR1

data set | $\rm ZCONF=1$ | $\rm ZCONF>1$ | Total
---|---|---|---
MOSAIC | 214 | 15.6% | 110 | 8.0% | 324 | 23.6%
UDF-10 | 35 | 19% | 17 | 9.2% | 52 | 28.3%
MXDF | 159 | 23.9% | 257 | 38.7% | 416 | 62.7%
All | 408 | 18.4% | 384 | 17.3% | 792 | 35.7%

Figure 25: Comparison of redshift statistics between the DR1 and DR2 data
releases. Light and dark colors are for low ($\rm ZCONF=1$) and high ($\rm
ZCONF>1$) redshift confidence, repectively.

####  6.1.3 Comparison with published spectroscopic redshifts

We also compared our redshifts with the compilation of all published
spectroscopic redshifts in the CDFS given by Merlin et al. (2021) prior to
2021242424We note that the MUSE DR1 redshifts published in Inami et al. (2017)
were not included in the Merlin et al. (2021) compilation.. The compilation
has 4951 entries, of which 292 are within the area of MUSE HUDF, which is 14
times smaller than the CDFS area.

A large majority (255 or 87%) of the published spectroscopic redshifts were
successfully matched to MUSE redshifts, with the exception of 26 sources (9%)
that were assigned a different redshift. Excluding the three galaxies with a
low MUSE redshift confidence (MID-7951, MID-8025 and MID-8038) for which the
proposed published solution is a viable alternative, we are confident that the
MUSE assigned redshift is the correct solution. In addition, 11 galaxies (4%)
of the published sample did not get any redshift in MUSE. We verify the
proposed redshift solution and confirm the absence of clear features in the
MUSE wavelength range. For these galaxies, the main feature (e.g., Ly emission
or absorption, [O iii]4959,5007 emission) lies outside the MUSE wavelength
range.

The scaled deviation of the published redshift with respect to MUSE, that is,
$\rm\Delta z=(z_{MUSE}-z_{Pub})/(1+z_{MUSE})$, has a standard deviation of
$\rm\Delta z=0.014$ after clipping 25 sources at $5\sigma$. We note that the
redshift errors reported in the literature are often greatly underestimated:
on average, they are 6 times smaller than the estimated error based on the
MUSE reference redshifts.

Just looking at the raw numbers (i.e., 2221 versus 292 spectroscopic
redshifts), one can appreciate the dramatic progress achieved by MUSE. As
shown in Fig. 26, the difference is even more striking at high redshift
($z>2.9$) where only 20 sources out of 1291 have previously been identified.
The same holds for faint magnitudes: most of the published spectroscopic
redshifts are from galaxies brighter than AB 24 in F775W, while the magnitudes
of MUSE sources extend to AB 30 and beyond (Fig. 21 left panel).

Figure 26: Comparison of the spectroscopic redshift distributions of ground-
based MOS and HST slitless GRISM published spectroscopic surveys in the HUDF
area and the MUSE DR2 redshifts.

####  6.1.4 Comparison with photometric redshifts

A comprehensive comparison of photometric redshifts with MUSE DR1 redshifts
was performed by Brinchmann et al. (2017). Photometric redshifts are essential
ingredients for on-going and upcoming cosmological surveys (e.g., KiDS,
Hildebrandt et al. 2017; LSST, Ivezić et al. 2019; Euclid, Laureijs et al.
2011; Roman, Spergel et al. 2015). Some surveys have very stringent
requirement in photo-z accuracy; for example, for weak lensing, the mean
redshift in each redshift bin must be constrained to better than $2\times
10^{-3}(1+z)$ (Newman et al., 2015). The goal of the Brinchmann et al. (2017)
study was to test the performance of photometric redshifts for galaxies in the
Hubble Ultra Deep Field down to the 30th magnitude.

One of the main conclusions of this study was the existence of systematic
biases at intermediate ($z=0.4-1.5$) and high redshift ($z>3$), and a
$\approx$10% proportion of outliers. We do not intend to repeat such a
detailed analysis here, but we are interested in whether these conclusions
based on DR1 are still valid for the more complete DR2 catalog.

For this comparison we use three different photometric redshifts catalogs. The
first catalog is from Rafelski et al. (2015) (hereafter BPZ-R15), which was
already used in Brinchmann et al. (2017). We use the Bayesian photometric
redshift values derived with the BPZ code (Benítez, 2000). The second catalog
contains the 3D-HST redshifts from Momcheva et al. (2016). These redshifts
were computed from a simultaneous fit of multiband HST photometry and 2D
spectra extracted from the HST WFC3 G141 GRISM. Although it is not a pure
photometric redshift catalog, we included it in the comparison given that the
additional GRISM data should provide improved redshifts with respect to purely
photometric redshifts. We note, however, this is only true in the wavelength
range covered by G141 and only for specific redshifts for which spectral
features (especially breaks) are sampled by the G141 wavelength range. It is
also restricted to object bright enough (i.e., with F140 JH magnitudes
brighter than 26). The last one is the recent Merlin et al. (2021) ASTRODEEP-
GS43 spectro-photometric catalog. It is based on HST CANDELS data complemented
with additional ground-based photometric data, reaching a total of 43
different photometric wide and narrow bands. Their photometric redshifts are
derived with three different codes (LePhare; Ilbert et al. 2006, EAzY; Brammer
et al. 2008, z-phot; Fontana et al. 2000), using the median of the three
derived values.

The input catalog is composed of 1711 sources with good and high confidence
redshifts ($\rm ZCONF>1$). The BPZ-R15, 3D-HST, and ASTRODEEP catalogs have
respectively 1404, 1088, and 1009 matched sources. As in Brinchmann et al.
(2017), we compute the normalized redshift error as

$\Delta z=\frac{z_{\rm MUSE}-z_{\rm phot}}{1+z_{\rm MUSE}}.$ (6)

We evaluated the evolution of $\Delta z$ with redshift and magnitude by
computing the median error and the outlier fraction in bins of about 100
elements. We note that the distributions of redshift and magnitude are
different for the three catalogs. For example, BPZ-R15 has more sources than
the other two catalogs, it extends to fainter magnitudes and has a larger
number of sources at high-z. Consequently, the redshift and magnitude bins
differ between catalogs. In each bin, the outlier fraction is computed as the
fraction of sources with $\Delta z>n\,\sigma_{MAD}$ where $\sigma_{MAD}$ is
the median absolute deviation around the median, adjusted to match a standard
deviation for a Gaussian distribution (Eq. 2 of Brinchmann et al. 2017). A
value of $n=3$ is used for simple outliers and $n=10$ for catastrophic
outliers.

The results are shown in Fig. 27. The BPZ-R15 curve (blue curve in top left
panel) displays the same behavior with redshift as observed in DR1 (Fig. 5
left panel of Brinchmann et al. 2017), with a systematic negative offset
($\rm\Delta z\approx-0.03$) at $\rm z\approx 1.2$ and a positive offset
($\rm\Delta z\approx+0.05$) at $\rm z\approx 3$. We note, however, that the
increased offset at $\rm z>4$ observed in DR1 is not observed in our larger
sample. The median standard deviation measured in BPZ-R15 is 0.04, only
slightly higher than the DR1 estimation. The evolution as a function of
magnitude (blue curve and shaded area in the top right panel) displays a
similar trend as for DR1: the standard deviation increases from 0.03 at $\rm
AB<25$ to 0.06 at $\rm AB>28$. The observed fraction of outliers is 9.4% and
2.8% for the catastrophic outliers. These figures are consistent with the DR1
measurements.

The comparison with the ASTRODEEP photometric redshifts is shown in green in
Fig. 27. The additional photometric bands have improved the accuracy of the
photometric redshifts. There are almost no detectable systematic offset with
either redshift or magnitude. The scatter is also much reduced compared to
BPZ-R15. It is less than 0.02 at AB 26 and then increases to 0.04 at AB 28.5.
We observe, however, that the average fraction of outliers increases, with 16%
and 9% for the $\rm 3\sigma$ and $\rm 10\sigma$ outliers, respectively. We
note that these numbers (scatter and outlier fraction) are higher than those
quoted by Merlin et al. (2021). Their underestimation probably results from
the limited spectroscopic sample used (see Sect. 6.1.3) in their tests.

Surprisingly, the comparison with the 3D-HST catalog (orange curve and shaded
area in Fig. 27) reveals a net degradation of performance at $\rm AB>26$. We
note that the Momcheva et al. (2016) catalog only includes grism information
for sources with JH F140 magnitudes brighter than 26. The increased scatter
must therefore be due to the photometric redshift estimate.

Figure 27: Redshift offﬀset between spectroscopic redshifts from MUSE and the
photometric redshifts from Rafelski et al. (2015) BPZ-R15 (blue), Momcheva et
al. (2016) 3D-HST (orange) and Merlin et al. (2021) ASTRODEEP (green)
catalogs, all normalized by $\rm 1+z_{MUSE}$. The top left panel show the
median (solid line) as a function of redshift and the scatter on the median is
shown by the shaded areas. The middle and bottom left panels show the fraction
of outliers at respectively $\rm 3\sigma$ and $\rm 10\sigma$. The right panels
display the same quantities as a function of HST F775W AB magnitude.

###  6.2 Sources without HST counterparts

During the inspection process (Sect. 5.7.4) we identified 424 (20%) MUSE
sources that cannot be matched to any of the four HST catalogs: R15 UVUDF
(Rafelski et al., 2015), CANDELS v2 (Whitaker et al., 2019), 3D-HST (Skelton
et al., 2014; Momcheva et al., 2016), ASTRODEEP (Merlin et al., 2021). In a
some cases (86), the source can still be clearly seen in the HST image. An
example of such ”missed” sources is shown in Fig. 28. This demonstrates the
limitations of source de-blending based solely on morphology and broadband
colors.

The vast majority of sources not found in the HST catalogs were classified as
”faint” or ”undetect” by the experts. In order to give more quantitative
ground for this classification, we performed our own HST photometry at the
locations of MUSE sources unmatched to the HST catalogs.

As in Maseda et al. (2018), we measure the S/N ratio for objects using the
HST/ACS mosaics from the XDF survey (Illingworth et al., 2013). We determine
the fluxes in apertures of 0.4 arcsec diameter centered on the DR2 catalog
positions. This aperture corresponds to a physical size of 3.1 kpc (2.2 kpc)
at $\rm z=2.9$ (6.6). The local background level is calculated by measuring
the standard deviation of the fluxes in 250 identical apertures spread
randomly within a $10\arcsec\times 10\arcsec$ cutout centered on the object,
with other objects masked according to the R15 segmentation maps. If the
aperture flux is greater than three times the local background level in an HST
band, then we consider the object “detected” in that band. Otherwise we use
the 3-$\sigma$ upper-limit to the flux in that band.

The computed magnitudes of these sources are shown in the right panel of Fig.
29. Using the $\rm 3\sigma$ S/N detection limit, we have 189 sources that fall
below the detection limit in all HST bands. We note that the majority of this
sample (90% or 380 sources) have AB magnitude greater than 29. In practice, as
shown in Sect. 6.1.4, no reliable photometric redshifts can be obtained for
these faint sources, even at the UDF HST depth (Fig. 27).

We observe that almost all sources without HST counterparts are
$\mathrm{Ly}\alpha$ emitters (left panel of Fig. 29). The nature of these
sources has been discussed by Maseda et al. (2018, 2020). The authors conclude
that they constitute the high equivalent-width tail of the $\mathrm{Ly}\alpha$
emitter population and that they have on average a high ionizing photon
production efficiency. Although their study is based on 103 faint sources only
and does not include the deepest part of the DR2 catalog related to MXDF
observations, these conclusions are not invalidated by the present data set.

Figure 28: Example of a MUSE source missed in HST catalogs. The upper left
panel display the HST F775W image ($2\aas@@fstack{\prime\prime}5\times
2\aas@@fstack{\prime\prime}5$) at the location of RID-21170. The red line
shows the R15 segmentation contour, while the yellow line displays the CANDELS
v2 segmentation contour of the corresponding source (ID 95547) in Whitaker et
al. (2019). Although all HST catalogs identify a single source, MUSE
unambiguously identifies two sources: MID-1319, an [O ii] emitter at $z=1.292$
(bottom central panel) and MID-6506 (bottom right panel), a
$\mathrm{Ly}\alpha$ emitter at $z=3.798$. The MUSE reference spectrum is shown
in the bottom left panel. From the corresponding [O ii] and Ly narrowband
images, one can see that the MID-1319 [O ii] emitter matches the brighter part
of the HST image, while the MID-6506 $\mathrm{Ly}\alpha$ emitter coincides
with the southeast secondary HST peak. Observed wavelength and flux units are
respectively Å and $\mathrm{10^{-20}\,erg\,s^{-1}\,cm^{-2}\,\AA^{-1}}$.
Figure 29: Statistics of MUSE sources without HST counterpart. Left panel:
comparison of the redshift distributions of MUSE sources with (brown color)
and without (blue color) HST counterparts. Right panel: estimated AB aperture
magnitude for MUSE sources without HST counterparts. The histogram is colored
according to the HST filters. Successful photometric measurements ($\rm
S/N>3$) are displayed in dark colors, while measurements with $\rm S/N<3$ are
displayed with light colors. In the latter case the reported magnitude is the
measured noise RMS standard deviation.

###  6.3 Completeness

The completeness measures how complete our sample is with respect to one of
its observational or physical properties (e.g., source magnitude or line
flux). It is an essential measurement to infer the LF of a given galaxy
population. To estimate the completeness, we have to evaluate the number of
missed sources at a given luminosity and therefore capture all the biases
introduced by the whole data chain (i.e., sample selection, observation, data
reduction, and analysis).

In this paper, the objective is to maximize the number of reliable redshifts
in the HUDF area. We therefore have used multiple detection methods: blind
detection with ORIGIN and HST prior extraction with ODHIN. We also used
additional available information (e.g., photometric redshifts) to choose
between various redshift solutions. Each detection method has its own bias.
While HST prior ODHIN detection is based on continuum detection and will
therefore be biased toward galaxies with stellar continuum, ORIGIN blind
detection will preferentially select high equivalent-width emission line
galaxies. Inferring the completeness of our catalog is therefore a challenging
task and beyond the scope of this paper.

A useful estimate of completeness can nevertheless be made by comparing our
catalog with the HST spectro-photometric catalogs. Given the great depth
reached by Hubble in the UDF (i.e., $\rm 5\sigma$ 29.5 AB F775W magnitude),
the HST catalogs should be quite complete and could be used as an unbiased
parent sample of broadband-selected galaxies. A rough estimate of HST Ly
completeness can be obtained for the $\mathrm{Ly}\alpha$ emitter population by
counting the fraction of MUSE sources with HST counterparts with respect to
the full $\mathrm{Ly}\alpha$ emitters sample. As shown in Fig. 30, the HST
catalogs can be considered as complete (90%) up to AB 29.5. At fainter
magnitudes, the HST catalogs have a significant drop.

We performed the comparison between MUSE and the R15 spectro-photometric
catalog in two areas: the MXDF region limited to the 100+ hours of depth and
the MOSAIC area with 10+ hours of depth, after removing the MXDF area. The
former covers an area of 0.85 arcmin2 and has 787 sources identified252525We
exclude all sources without measurable F775W magnitude. in the R15 catalog,
the latter covers 7.48 arcmin2 and has 3619 sources identified in R15. The
average depth achieved with MUSE in these two areas is 140.8 and 10.8 hours,
respectively. In each sample, we measure the fraction of HST sources that have
a MUSE redshift in a given magnitude or redshift interval. The completeness
evolution with respect to magnitude and photometric redshift is given in Fig.
31.

For a given magnitude, the evolution of completeness with redshift is not
flat: first, we see a clear drop in the $z=2-3$ range, corresponding to the
redshift range of the MUSE desert, then we measure a higher completeness at
$z>3$, indicating that the fraction of $\mathrm{Ly}\alpha$ emitters is high in
this redshift range. In the MOSAIC sample at a 10-hour depth, the 50%
completeness is achieved at a magnitude of 26.6 (F775W) in the redshift range
$z=0.8-1.6$ and 27.5 at $z\approx 4$. In the deepest MXDF observations, the
corresponding magnitude for the 50% completeness are 27.6 and 28.7 for
$z=0.8-1.6$ and $z=3.2-4.5$, respectively.

It is instructive to compare these figures with those obtained for the
published spectroscopic redshifts sample (Merlin et al. 2021 compilation). The
50% completeness rate is achieved at a magnitude of 25.2 (F775W) in the $\rm
z<3$ redshift range. At higher redshifts ($\rm z>3$) the completeness drops
very rapidly with increasing magnitude: 20% for AB 26, 5% for AB 27.

Figure 30:  Fraction of MUSE $\mathrm{Ly}\alpha$ emitters with HST
counterparts as a function of AB F775W magnitude in bins of 0.5 magnitude.
Figure 31:  DR2 completeness with respect to HST deep imaging sample.
Completeness is computed as the number of galaxies with a MUSE redshift in a
0.31 HST AB F775W magnitude (left panel) or 0.5 photometric redshift interval
(right panel). The measurements have been performed for the deep MXDF
(+100-hour depth, orange bars) and MOSAIC (+10-hour depth, blue bars) data
sets and corresponding areas.

###  6.4 Physical properties

In this section we derive the stellar masses and star-formation rates (SFRs)
for our galaxy sample using spectral energy distribution (SED) fitting. Since
the wavelength range of MUSE does not extend far enough into the infrared to
provide sufficient constraints for the SED fit, we selected the subsample of
1664 galaxies (75%) with $\rm ZCONF>0$ that have been identified in the HST
R15 catalog. The R15 catalog has 11 photometric bands262626WFC3/UVIS F225W,
F275W, and F336W; ACS/WFC F435W, F606W, F775W, and F850LP and WFC3/IR F105W,
F125W, F140W, and F160W. ranging from the NUV ($0.21\mu m$) to the WFC3/IR
($1.5\mu m$).

We used the high-z extension of the code Magphys (da Cunha et al., 2008, 2015)
with a minimum stellar mass of $\rm 10^{6}M_{\odot}$ (see Sect. 3.2 of Maseda
et al. 2017). The redshift of each object is fixed to the MUSE spectroscopic
redshift. In cases where the fitting results in un-physically small errors in
stellar mass and/or SFR (i.e., when the models cannot accurately reproduce the
photometric data points), we adopted a fixed error of $\pm 0.3$ dex (the
typical uncertainty in determining stellar masses from fits to broadband
photometry; e.g., Conroy 2013; Maseda et al. 2014).

In addition to Magphys, we further utilized the Prospector SED fitting code
(Johnson et al., 2021) to derive stellar masses and SFRs for the same R15
photometric catalog. Through Prospector we are able to consider the
contribution from emission lines to galaxy photometry. Given that a large
fraction of our sources have strong emission lines, constraining the
contribution of the emission lines in deriving galaxy properties is important.

In Prospector we used FSPS models (Conroy et al., 2009) with the MILES
spectral library and MIST isochrones to generate synthetic spectra to match
with the observations. These models include contributions from nebular
continuum and emission as described by Byler et al. (2017). We used a
nonparametric approach to fit the observed photometry using the
“continuity_sfh” template in Prospector. This computes the stellar mass
produced in fixed time bins with a smoothness prior (Leja et al., 2019). We
defined 7 time bins to fit for stellar mass. The first two bins were fixed at
0-30 and 30-100 Myr in lookback time. The furthest bin was fixed to be between
85% - 100% in lookback time defined by the MUSE spectroscopic redshift. The
remaining four bins were evenly distributed between the second and final bin
in logarithmic space. We used a Calzetti et al. (2000) dust attenuation law
and allowed the V-Band dust optical depth to vary as a free parameter between
0 and 2.0. Gas-phase metallicity and stellar metallicity were tied together
and were allowed to vary between 0.01 and $\approx$1.5 times the solar
metallicity. We let the gas ionization parameter vary between $\rm log(U)=-4$
and $-1$. In total we had 10 free parameters (6 SFH bins + stellar mass, dust,
metallicity, and ionization parameter). Parameters were sampled using the
Dynamic Nested Sampling method presented by Speagle (2020) using 200 live
points and batches and 50 nested walks.

The stellar masses and the SFR averaged over the last 100 Myr are given in the
main catalog (Sect. A.2). The error for the masses and SFR is based on the
1-sigma percentiles of the nested sampling chains. The stellar mass output by
Prospector is the “total formed mass.” In order to convert this to the
observed mass, we use the mass correction factor from the maximum a posteriori
solution (the posterior sample with the highest posterior possibility) to
derive an approximate correction factor for the marginalized stellar mass
values.

Two examples of the SED fit are given in Fig. 32. In the case of MID-3 (upper
left panel), it can be seen that the SED fit of Magphys, unlike Prospector,
overestimates the continuum due to the presence of strong emission lines272727
There are also some small differences in the photometry between MUSE and HST,
as discussed in Sect. 4.1, that may account for a small offset between the two
sets of measurements.. In Fig. 33, we compare the Magphys and Prospector
stellar mass and SFR results for the full sample. We note that Prospector
tends to derive higher masses than Magphys with a median offset of 0.25 dex.
This is a known feature of Prospector that is documented in Leja et al. (2019,
2020). The stellar mass and SFR distributions are presented in the first two
panels of Fig. 34. As expected, the MUSE catalog probes low mass star forming
galaxies with a median mass of $\rm 6.2\times 10^{8}\,M_{\odot}$ and a median
SFR of $\rm 0.4\,M_{\odot}yr^{-1}$. The last panel of Fig. 34 displays the
traditional star forming main sequence. We note that our main sequence
deviates significantly from the best-fit star forming main sequence at $z=1$
derived from a compilation of 25 studies by Speagle et al. (2014) (dashed blue
line). Such an offset between Prospector SED results and previous studies has
been identified and discussed by Leja et al. (2021). As noted by the authors
there are many reasons for the difference in inferred stellar mass and SFR
between the Prospector’s nonparametric star-formation histories and the
classical approaches (see their Sect. 6.2 for a detailed discussion of the
origin of this difference). We also compared our results at $z=1$ with those
of Boogaard et al. (2018) (Fig. 34 solid blue line). The Boogaard et al.
(2018) study was based on DR1, still using HST photometry SED to estimate
stellar mass but with a different derivation of SFR based on the H and H
recombination line flux. As shown in the plot, there is still a small offset
between our estimate and the results from Boogaard et al. (2018), but it is
much smaller compared to the values of Speagle et al. (2014). This residual
offset is most likely due to the different estimate of the SFR. Further
investigations will be needed, for example into the completeness and SED
degeneracy, to interpret this main sequence, but such an analysis is beyond
the scope of this paper.

Finally, we recall that this SED analysis does not take into account the
population of 577 (25%) galaxies without an HST counterpart in the R15
catalog. The majority of these sources are high equivalent-width
$\mathrm{Ly}\alpha$ emitters that are too faint in the HST broadband deep
images to get meaningful SED measurements. However, as demonstrated by Maseda
et al. (2018, 2020), this selection yields on average faint
($\mathrm{M_{UV}\approx-15}$), blue ($\mathrm{\beta\approx-2.5}$) star-forming
galaxies. Although we cannot directly determine the stellar masses from the UV
continuum alone, an extrapolation to the Duncan et al. (2014)
$\mathrm{M_{UV}-M_{\odot}}$ relation to $\mathrm{M_{UV}}$ = $-$15 implies that
these galaxies should have stellar masses below $\mathrm{10^{7}M_{\odot}}$ at
these redshifts.

Figure 32: Two example SED fits: the $z=0.62$ MID-3 [O ii] emitter (top
panels) and the $z=4.36$ MID-1283 $\mathrm{Ly}\alpha$ emitter (bottom panels).
The left column shows the MUSE spectrum (in black) and the Prospector (blue
line) and Magphys (orange line) SED fits. The right column displays the SED
fits and the HST R15 photometry (black symbols). The shaded area in the right
column indicates the MUSE wavelength range. Flux and wavelength units are
$\mathrm{10^{-20}\,erg\,s^{-1}\,cm^{-2}\,\AA^{-1}}$ and $\rm\mu m$,
respectively.  Figure 33: Comparison of the Prospector and Magphys stellar
mass and SFR results for the SED fits of all MUSE galaxies identified in the
R15 catalog.  Figure 34: Stellar mass (top left panel) and SFR (top right
panel) distributions for the MUSE galaxies identified in the R15 catalog.
Values are given for the Prospector SED fit in the $\approx$100 Myr lookback
time. The stellar mass and the SFR are in $\rm log\,M_{\odot}$ and $\rm
log\,M_{\odot}\,yr^{-1}$ units, respectively. The vertical black lines show
the median values of the distributions. The bottom panel displays the star
forming galaxies main sequence $\rm log\,SFR(log\,M)$. Symbols are colored
according to the redshift. The z=1 main sequence derived by Boogaard et al.
(2018) and Speagle et al. (2014) are shown as blue solid and dashed lines,
respectively. See the main text for a discussion of the offset from our
inferred values.

###  6.5 Peculiar objects

Among this rich set of data, we give in what follows some examples of
interesting objects. The first category is the $\mathrm{Ly}\alpha$ emitters
with very extended Ly emission. These objects are often called Lyman alpha
blobs (LABs) in the literature. They are rare and mainly found in overdense
regions like the two LABs discovered in the SSA22 supercluster (Steidel et
al., 2000). While most $\mathrm{Ly}\alpha$ emitters and rest-frame UV-selected
galaxies exhibit diffuse Ly emission (Wisotzki et al., 2016; Leclercq et al.,
2017; Kusakabe et al., 2022), LABs differ in terms of their high luminosity
(typically brighter than $\mathrm{10^{43}erg\,s^{-1}}$) and their large size
(typically larger than 50 kpc; e.g., Herenz et al. 2020). A recent review is
provided in Kimock et al. (2021). We show in Fig. 35, three examples of
$\mathrm{Ly}\alpha$ emitters with such extended Ly emission. The first object
(MID-208 and MID-6329) is a $z=3.32$ source detected in the UDF-10 at a
29-hour depth. With a Ly luminosity of $\mathrm{{8}\times
10^{42}\,erg\,s^{-1}}$ and a size of $\approx$60 kpc, it enters the LAB
category. This object was already reported by Vanzella et al. (2017) using 20
hours of the same GTO observations. The second example (MID-1056 and MID-3621)
was discovered in the MOSAIC at an 11-hour depth. This $z=3.06$ LAB is also
bright ($\mathrm{{1.4}\times 10^{43}\,erg\,s^{-1}}$ Ly luminosity) and
extended ($\approx$60 kpc). This object coincides with a luminous type II
active galactic nucleus identified in the Chandra 7Ms catalog (Luo et al.,
2017). The extended Ly emission was already reported by den Brok et al. (2020)
using the same data set. It is part of a 1.2 Mpc long Ly filament (Bacon et
al., 2021). The third example, the $z=3.2$ MID-1530 $\mathrm{Ly}\alpha$
emitter, is less extended (30 kpc) and luminous ($\mathrm{{1.4}\times
10^{42}\,erg\,s^{-1}}$ Ly luminosity) than the other two and would not be
formally categorized as a LAB. It is representative of the class of the large
Ly halos found around $\mathrm{Ly}\alpha$ emitters (e.g., Leclercq et al.
2020).

In addition to the MID-1056/3621 type II active galactic nucleus, we have two
type I QSOs in the sample (Fig. 36). The first is MID-977, a $z=1.2$ object
with broad Mg ii emission and some absorbers detected along the line of sight.
The second is MID-1051, a $z=3.2$ QSO with damped Ly absorption and C iv, He
ii, and C iii] emission.

As previously stated, $\mathrm{Ly}\alpha$ emitters are typically surrounded by
a Ly halo that is detected out to approximatively ten times the continuum size
(e.g., Leclercq et al., 2017). Moreover, the peak of the Ly emission can be
offset with respect to the UV continuum peak (e.g., Ribeiro et al., 2020;
Claeyssens et al., 2022). We present an extreme example of this in Fig. 37,
where the $z=4.8$ $\mathrm{Ly}\alpha$ emitter (MID-1264) has a 4.1 kpc (0.6
arcsec) offset between the UV and Ly narrowband centers. We note that such an
offset is much larger than the $0.6\pm 0.05$ kpc average offset observed by
Ribeiro et al. (2020) in their deep VIMOS slit data. Nevertheless, these
offset measurements rely on a correct match between the UV continuum source
(in our case the deep HST image) and the $\mathrm{Ly}\alpha$ emitter. This is
not always an easy task, and MID-1264 is a good example. The
$\mathrm{Ly}\alpha$ emitter was matched to RID-6956. This object is bright
enough (25.6 AB F775W magnitude) to have a reliable photometric redshift. All
HST catalogs give similar estimates, with for example $z_{\rm phot}=4.77\pm
0.07$ for the ASTRODEEP catalog. Our measurement of $z=4.78$ fits well within
the photo-z error bar, giving confidence that the $\mathrm{Ly}\alpha$ emitter
is matched with the HST source. In addition, all other nearby HST sources are
much farther from the location of the Ly emission peak and their MUSE or
photo-z redshifts indicate that they are foreground objects. However, the high
S/N of the Ly line profile shows evidence for a small dip at 7029 Å that can
be interpreted as the superposition of two $\mathrm{Ly}\alpha$ emitters at
$z=4.7806$ and $z=4.7885$ respectively. In this case, the fainter Ly line
would be associated with RID-6956, while the brighter Ly line would be a
$\mathrm{Ly}\alpha$ emitter without an HST counterparts, satellite to the main
galaxy.

There are many examples where the HST sources have no reliable photo-z, are
located in a crowded environment or have a complex morphology. In these cases,
matching Ly emission to HST imaging is a challenge. In general, we tend to
favor HST-undetected $\mathrm{Ly}\alpha$ emitter when the offset is large
enough. This is in agreement with a recent study that used a large sample of
lensed $\mathrm{Ly}\alpha$ emitters observed with MUSE (Claeyssens et al.,