# An optimal ALMA image of the Hubble Ultra Deep Field in the era of JWST:
obscured star formation and the cosmic far-infrared background

Ryley Hill,1 Douglas Scott,1 Derek J. McLeod,2 Ross J. McLure,2 Scott C.
Chapman3,4,1 and James S. Dunlop2
1Department of Physics and Astronomy, University of British Columbia, 6225
Agricultural Road, Vancouver, BC, V6T 1Z1, Canada
2Institute for Astronomy, University of Edinburgh, Royal Observatory,
Edinburgh, EH9 3HJ, UK
3Department of Physics and Atmospheric Science, Dalhousie University, Halifax,
NS, B3H 4R2, Canada
4NRC – Herzberg Astronomy and Astrophysics, 5071 W Saanich Rd, Victoria, BC,
V9E2E7, Canada E-mail<EMAIL_ADDRESS>
(September 2023)

###### Abstract

We combine archival ALMA data targeting the Hubble Ultra Deep Field (HUDF) to
produce the deepest currently attainable 1-mm maps of this key, deep,
extragalactic survey field. Combining all existing data in Band 6, our deepest
map covers $4.2\,{\rm arcmin}^{2}$, with a beamsize of $1.49\,{\rm
arcsec}\,{\times}\,1.07\,{\rm arcsec}$ at an effective frequency of $243\,$GHz
or $\lambda=1.23\,$mm. It reaches a minimum pixel rms of $4.6\,\mu{\rm
Jy}\,{\rm beam}^{-1}$, with an area of 1.5 arcmin2 reaching below 9.0 $\mu$Jy
beam-1. This is an improvement of at least 5 per cent compared with the best
previously published map over the same area, with as much as a 50 per cent
improvement in some regions. We also make a wider, but shallower map, covering
$25.4\,{\rm arcmin}^{2}$. We detect 45 galaxies in the deep map down to
$3.6\,\sigma$, including 10 more 1-mm sources than previously detected within
the same area. 38/45 of these galaxies have a JWST identification from the
JADES NIRCam imaging and we find that the new sources are typically faint and
red. A stacking analysis on the positions of ALMA-undetected JADES galaxies as
a function of stellar mass and redshift yields significant detections in our
image for redshifts below 4 and stellar masses from $M_{\ast}=10^{10.4}\,{\rm
M}_{\odot}$ down to $M_{\ast}=10^{8.4}\,{\rm M}_{\odot}$, enabling us to
extract a modest amount (${\simeq}\,10\,$per cent) of additional stacked
signal from our map as compared to previous analyses. Combining detected
sources and statistical estimates from stacking, we have detected
$(10.0\,{\pm}\,0.5)\,{\rm Jy}\,{\rm deg}^{-2}$ of the cosmic infrared
background (CIB) at 1.23 mm. Although this is short of the (uncertain) average
background level of about 20 Jy deg-2, after taking into account intrinsic
fluctuations in the CIB, we find that our measurement is consistent with the
average background level if the HUDF is a mild (${\simeq}\,2\,\sigma$)
negative fluctuation. This suggests that within the HUDF, JWST may have
detected essentially all of the galaxies that contribute to the CIB at 1.23
mm, although we have not yet directly detected all of these galaxies with
ALMA. From our stacking analysis we predict that the field contains around 60
additional galaxies with 1.23 mm flux densities averaging around 15 $\mu$Jy,
and over 300 galaxies with flux densities of a few $\mu$Jy. However, the
contribution of these fainter more modestly-obscured objects to the background
is small, and converging, as anticipated from the now well-established strong
correlation between galaxy stellar mass and obscured star formation.

###### keywords:

methods: data analysis – techniques: interferometric – galaxies: formation –
galaxies: starburst – submillimetre: galaxies

††pubyear: 2023††pagerange: An optimal ALMA image of the Hubble Ultra Deep
Field in the era of JWST: obscured star formation and the cosmic far-infrared
background–15

##  1 Introduction

The Hubble Ultra Deep Field (HUDF) is probably the most well-studied
extragalactic region of the sky, containing some of the deepest optical and
near-IR exposures obtained to-date (e.g. Beckwith et al., 2006; Illingworth et
al., 2013; Koekemoer et al., 2013; Eisenstein et al., 2023). Studies of the
HUDF can tell us about how galaxies have evolved from the earliest times to
the present day. We now know that a significant fraction of the cosmic star
formation occurred within distant galaxies full of dust (e.g. Casey et al.,
2015; Koprowski et al., 2017), making them effectively invisible at optical
and near-IR wavelengths, for all but the deepest images. Yet these galaxies
are bright at millimetre (mm) and submillimetre (submm) wavelengths, and so in
order to complete our understanding of galaxy evolution, we must also survey
the HUDF in this waveband.

The best telescope at these wavelengths available today is the Atacama Large
Millimeter/submillimeter Array (ALMA). ALMA has been used several times to
survey the HUDF at wavelengths around 1 mm (Dunlop et al. 2017; ASAGAO,
Hatsukade et al. 2018; ASPECS, González-López et al. 2020; GOODS-ALMA, Gómez-
Guijarro et al. 2022), and it has also been used to follow up interesting
individual targets within the HUDF (e.g. Fujimoto et al., 2017; Cowie et al.,
2018). Yet in contrast to the tens of thousands of galaxies detected in the
optical, only a few dozen mm-bright galaxies have been found in these follow-
up observations so far. Although often referred to as ‘SMGs’ (for
‘submillimetre galaxies’), since we are working here at wavelengths slightly
longer than 1 mm, we will refer to them with the more generic acronym DSFG
(for ‘dusty star-forming galaxy’).

In order to detect more galaxies around 1 mm with ALMA, we can combine all of
the data into a single image. In a similar way, two decades ago, a ‘super-map’
of the GOODS-N field was constructed by combining SCUBA data sets at 850
$\mu$m (Borys et al., 2003), and multiwavelength counterparts from the Hubble
Space Telescope (HST) were identified and used to study the properties of the
detected submm-luminous sources (Pope et al., 2005). With this same motivation
we have undertaken an archival project to combine all of the ALMA observations
of the HUDF taken around 1 mm, with the goal of finding new DSFGs, identifying
their counterparts, assessing how much of the background light has been
resolved and providing an image that can be used by others for stacking
analyses.

In Section 2 we describe how we retrieve the data from the ALMA archive and
combine it into single continuum images in $uv$ space. In Section 3 we outline
an alternative approach to combining data subimages in real space, which
provides a test of the reliability of our $uv$-space combination. In Section 4
we provide our new galaxy catalogue and compare it to previous catalogues. In
Section 5 we describe how our results lead to a new estimate of the resolved
fraction of the cosmic infrared background (CIB) at 1 mm. We discuss in
Section 6 improvements in our data products and results compared to what was
previously available at 1 mm in the HUDF and we conclude in Section 7.
Appendix A presents a wider (and shallower) map and gives a supplementary list
of sources, while Appendix B provides cutouts of our ALMA sources overlaid
over JWST F356W imaging.

##  2 Data retrieval and processing

###  2.1 Obtaining archival ALMA data

We focus on ALMA Band 6, which spans a wavelength (frequency) range of 1.1–1.4
mm (210–270 GHz). This band currently contains the most extensive ALMA
observations of the HUDF, and so this is where we expect to produce the
deepest archival map. To begin, we queried the ALMA
archive111https://almascience.nrao.edu/asax for all Band-6 observations
centred at 03h32m39.0s $-$27∘47′29.1′′ and overlapping within a radius of 1.5
arcmin. There are a total of 12 unique programmes that satisfy this criterion,
and we selected all of these for our combined map. For some programmes, only a
fraction of the total time was spent on the HUDF (e.g. for individual
pointings of targeted galaxies or for large surveys extending beyond the
HUDF), and for these cases we only selected the data that overlap with our
region of interest. For more details see Table 1, where we summarize all of
the data used here to produce the combined map. We also downloaded the raw
data from the ALMA archive centred at the same position, but extending out to
3.5 arcmin, which we used to make a shallower but larger map; there are a
total of seven additional programmes satisfying this criterion, and these are
listed in Table 4.

For each of the observations, we downloaded the raw $uv$ data and calibrated
it using the provided ScriptForPI and the CASA (McMullin et al., 2007) version
used by the observatory at the time of the observations. We split the science
targets from the calibration targets, then time-averaged the data by 30 s and
averaged the frequency channels by a factor of 4 in order to reduce the volume
of data. These tasks were carried out using the Canadian Advanced Network for
Astronomy Research (CANFAR) platform (Major et al., 2019), which provides easy
access to all versions of CASA, as well as ample storage space and memory.

###  2.2 Obtaining archival ALMA images

In addition to downloading the $uv$ data for each programme, we also obtained
the imaging products made available by the observatory. For every tuning of
every target given in Table 1, the observatory provides a single image made
using the multi-frequency synthesis (MFS) mode, where visibilities in each
channel are mapped to a single $uv$-plane, and therefore represent the mean
value of the sky at a characteristic frequency (usually the central frequency
of the channels) weighted by a spectral function. For programmes carried out
in earlier ALMA cycles, we use the Additional Representative Images for Legacy
(ARI-L; Massardi et al., 2021), which are reduced in a way more similar to
later cycles. These images are primary-beam-corrected, but the primary-beam
image is also available for download. The final number of HUDF images
downloaded from the ALMA archive is 61, along with their corresponding 61
primary-beam images.

Table 1: ALMA projects downloaded from the ALMA archive and used to create the
final combined 1-mm image.

Project ID | Target name(s) | Frequency range | Map rmsa | Synthesized beamb | Sky coveragec
---|---|---|---|---|---
|  | [GHz] | [$\mu$Jy beam-1] | [arcsec$\,{\times}\,$arcsec] | [arcmin2]
2012.1.00173.Sd | HUDF | 211.2–231.2 | 24 | 0.62$\,{\times}\,$0.52 | 7.05
2013.1.00718.Se | UDF1 | 212.2–272.0 | 15 | 1.77$\,{\times}\,$0.90 | 1.15
2013.1.01271.S | UDF6462 | 223.3–244.1 | 9 | 0.33$\,{\times}\,$0.26 | 0.31
2015.1.00098.Sd | HUDF-JVLA-ALMA | 244.3–271.8 | 57 | 0.20$\,{\times}\,$0.16 | 34.7
2015.1.00543.Sd | GOODS-S | 255.1–274.7 | 130 | 0.26$\,{\times}\,$0.22 | 54.2
2015.1.00664.S | KMOS3DGS4-24110 | 273.2–274.7 | 91 | 0.15$\,{\times}\,$0.14 | 0.22
2015.1.01096.S | UDF-640-1417 | 222.0–255.2 | 11 | 0.99$\,{\times}\,$0.74 | 0.29
2015.1.01447.S | UDF0 | 212.0–272.0 | 6 | 1.24$\,{\times}\,$0.91 | 0.29
2015.A.00009.S | Bo15Hz27 | 226.4–245.0 | 13 | 0.73$\,{\times}\,$0.46 | 0.30
2016.1.00324.Ld | UDF_mosaic_1mm | 212.1–272.0 | 12 | 1.37$\,{\times}\,$0.93 | 3.72
2017.1.00755.Sd | GOODS-S | 255.1–274.9 | 120 | 1.31$\,{\times}\,$0.82 | 54.1
2018.1.00567.S | ASAGAO27, 35, 40, 45 | 244.3–262.9 | 24 | 0.65$\,{\times}\,$0.47 | 1.03

  * •

a Approximate map rms estimated using the observatory MFS data products. For
programmes with a single tuning, this is the rms of the primary-beam-
uncorrected map after masking all known sources. For programmes with multiple
tunings, we follow the same procedure and then estimate the rms of the
weighted mean of the images as $1/\big{(}\Sigma_{i}1/\sigma_{i}^{2}\big{)}$.
The actual rms from combining images at different tunings in $uv$ space using
the MFS mode will generally be less than this estimate.

  * •

b Maximum synthesized beam across all frequencies of a given data set. The
actual synthesized beam from combining data taken in different tunings using
the MFS mode may differ slightly from this estimate.

  * •

c The sky coverage overlapping with our data selection criteria, which may be
less than the total sky coverage of the given programme.

  * •

d Programmes used to make a mosaic with uniform noise.

  * •

e The imaging products from this data set were not used to create the combined
image in the image plane because the beamsize is much larger than for the
other data sets.

###  2.3 Data combination in the uv plane

Our main goal is to produce the deepest possible map of the central HUDF
region, which essentially corresponds to the footprint of ASPECS (see
González-López et al., 2020). We do so by combining all of the available data
within this region. Our secondary goal is to produce the deepest possible map
of the wider HUDF region, corresponding essentially to the ASAGAO footprint
(see Hatsukade et al., 2018). To do this, we carried out the following
procedure.

All of the time-averaged and frequency-averaged $uv$ data were imaged using
the standard CASA task tclean. For the data presented in Table 1, the pixel
size was set to 0.2 arcsec. We ran tclean in MFS mode, which, as discussed
above, scales all of the $uv$ points to the average frequency of the data
being combined; here, the average frequency is about 243 GHz (or 1230 $\mu$m).
We also chose natural weighting (each $uv$ point is weighted by its
instrumental noise), but with a 250 k$\lambda$ $uv$ taper. We set the cutoff
of the map to 0.2 times the primary beam, which is where the footprint of our
map roughly matches the footprint of the ASPECS map. We also produced a map
with no $uv$ taper, setting the pixel size to 0.18 arcsec (because the
synthesized beam is somewhat smaller). Lastly, we cleaned the image down to 20
$\mu$Jy beam-1, similar to the cleaning level chosen by ASPECS. The final
synthesized beamsize of the $uv$-tapered map is $1.49\,{\rm
arcsec}\,{\times}\,1.07\,{\rm arcsec}$, while for the map with no $uv$ taper
the synthesized beamsize is $1.32\,{\rm arcsec}\,{\times}\,0.92\,{\rm
arcsec}$. The pixel rms in the tapered image (which we use for further
analysis in this paper) has a minimum value of 4.6 $\mu$Jy beam-1 in the
deepest region, while the rms is less than 9.0 $\mu$Jy beam-1 over an area of
1.5 arcmin2.

For the shallower, larger map, we exclude all of the observations targeting
the deep central region (this includes all of the ASPECS data, the data from
Dunlop et al. 2017, and two deep pointings towards the eastern corner of the
HUDF); this is to ensure that the synthesized beam of the larger map is not
dominated by data from the central region (see Table 4 for details). This map
is therefore mostly a combination of the ASAGAO survey data and the GOODS-ALMA
survey data (both high and low resolution). However, we include all of the
ASAGAO follow-up observations even if they lie within the deep central region
(programme ID 2018.1.00567.S), since they roughly match the array
configurations and frequencies of the ASAGAO survey and the GOODS-ALMA survey.
We set the pixel size to 0.06 arcsec for the map with no $uv$ tapering, and
0.12 arcsec for the map with 250 k$\lambda$ $uv$ tapering, and set the cutoff
of the map to 0.4 times the primary beam, which produced a map with a
footprint roughly similar to the ASAGAO footprint. All other tclean parameters
were kept the same as described above, except for the cleaning threshold,
which was increased to 100 $\mu$Jy beam-1 (similar to the level chosen by
ASAGAO). The final synthesized beamsize of the $uv$-tapered map is $0.87\,{\rm
arcsec}\,{\times}\,0.64\,{\rm arcsec}$, while for the map with no $uv$ taper
the synthesized beamsize is $0.36\,{\rm arcsec}\,{\times}\,0.27\,{\rm
arcsec}$. The pixel rms in the tapered image (which we also use for further
analysis in this paper) has a typical value of 60 $\mu$Jy beam-1 outside of
individual pointings, and 20 $\mu$Jy beam-1 within individual pointings.

Lastly, we produced maps with uniform noise properties by excluding programmes
with individual pointings targeting known galaxies, and instead focusing on
large survey programmes. Specifically, for our deep central mosaic we combined
data from the programmes 2012.1.00173.S, 2015.1.00098.S, 2015.1.00543.S,
2016.1.00324.L and 2017.1.00755.S, while for the shallower, larger map we only
used data from the programmes 2015.1.00098.S, 2015.1.00543.S and
2017.1.00755.S (see Tables 1 and 4, respectively). We used the same tclean
parameters as above, generating versions with no $uv$ tapering and with 250
k$\lambda$ $uv$ tapering. The final synthesized beamsizes were effectively
unchanged compared to the maps that included individual pointings.

For the remaining analysis in this paper we focus on the tapered maps with all
of the available data combined. The additional maps, with no tapering and
uniform noise properties, were used to check the reliability of our flux
densities, and we found no systematic differences between the measurements.

To calculate the corresponding noise maps for the combined images, we first
created a mask using catalogues of previously-detected galaxies. In particular
we took the catalogues from Dunlop et al. (2017), ASAGAO (Hatsukade et al.,
2018), ASPECS (González-López et al., 2020) and GOODS-ALMA (Gómez-Guijarro et
al., 2022). For each galaxy we masked a region 3 times larger than the
beamsize, and then multiplied this by the signal map. We next made cutouts of
5 times the beamsize around each pixel and calculated the rms. The resulting
noise map was then smoothed with a Gaussian kernel with the same size as the
cutout regions to remove artefacts.

All of these data products are made publicly available, including the primary-
beam-corrected maps, the noise maps, the primary beam maps and the synthesized
beam maps.222https://doi.org/10.5683/SP3/YWBVWH The final signal map, noise
map and signal-to-noise ratio (S/N) map for the deep central mosaic with a 250
k$\lambda$ taper, including individual pointings, is shown in Fig. 1, while
the same maps for the larger, shallower mosaic are shown in Fig. 13. The deep
central mosaic covers an area of 4.2 arcmin2 out to our primary beam threshold
of 0.2, while the large shallow mosaic covers an area of 25.4 arcmin2 out to
our primary beam threshold of 0.4.

###  2.4 Optical and infrared galaxy catalogues

The HUDF has been the target of extensive multiwavelength observations that we
can use to complement our 1.23-mm image. The Cosmic Assembly Near-IR Deep
Extragalactic Legacy Survey (CANDELS) catalogue of galaxies in the GOODS-S
field333https://archive.stsci.edu/prepds/candels (Guo et al., 2013) contains a
summary of much of the deep optical-to-near-IR imaging in the HUDF, and it is
expected that most of the ALMA-detected galaxies should have a counterpart in
this catalogue. Detections in the CANDELS catalogue were obtained from HST’s
WFC3 instrument in the F160W band, so it is a 1.5-$\mu$m-selected catalogue.
The entire catalogue covers a much larger region than our deep map, but the
deepest region of the catalogue covers roughly the same area.

The first public data release from the the JWST Advanced Deep Extragalactic
Survey (JADES) NIRCam survey of the HUDF is also now
available444https://archive.stsci.edu/hlsp/jades (Rieke & the JADES
Collaboration, 2023), in addition to the JWST Extragalactic Medium-band Survey
(JEMS, Williams et al., 2023). While the JADES survey is still ongoing and
will eventually be deeper, it already contains considerably more galaxies per
unit area than the CANDELS catalogue. JADES images range from 0.9 $\mu$m to
4.4 $\mu$m, and we chose to select our catalogue at $3.5\,\mu$m, as a
compromise between the longest possible wavelength and the depth. Our entire
deep map is covered by the JADES survey at approximately uniform sensitivity.
Throughout this paper we primarily make use of the JADES imaging and
photometry, with the CANDELS catalogue used to demonstrate new improvements
thanks to JWST.

####  2.4.1 JADES photometric redshifts and stellar masses

In order to construct photometry catalogues, we ran SourceExtractor (Bertin &
Arnouts, 1996) in dual-image mode, with NIRCam F356W used as the detection
image. All imaging was homogenized to match the point-source function (PSF) of
the F444W imaging in order to minimise any colour systematics arising from
differences in the PSF. We performed isophotal photometry on all available
JADES and JEMS NIRCam bands, as well as ancillary HST ACS imaging from the
Hubble Legacy Fields GOODS-South data release (see Illingworth et al., 2013;
Whitaker et al., 2019) and ground-based VIMOS $U$-band imaging over GOODS-
South (Nonino et al., 2009). To extract robust photometry from this lower-
resolution imaging, we utilised TPHOT (Merlin et al., 2015), using positional
and surface-brightness information from the higher-resolution NIRCam imaging
as input.

To calculate photometric redshifts and hence stellar masses, we performed
spectral-energy-distribution (SED) fitting using LePhare (Arnouts & Ilbert,
2011), with a Bruzual-Charlot (Bruzual & Charlot, 2003) template set,
incorporating a Chabrier (Chabrier, 2003) initial mass function. The template
set includes a Calzetti (Calzetti et al., 2000) dust-attenuation law, with
range of reddening $A_{V}=[0.0,6.0,0.2]$. We used two metallicities, 0.2 Z⊙
and Z⊙, and included exponentially declining star-formation histories with
$\tau$ ranging between 0.1 and 15 Gyr. Where the object has a known
spectroscopic redshift, we fixed to this redshift in obtaining the stellar
mass.

The median uncertainty in our photometric redshifts for all of our HUDF
galaxies is $\Delta z/(1\,{+}\,z)\,{\simeq}\,0.06$. For stellar masses, in
this paper we only care about relative uncertainties, since we only use this
parameter to sort galaxies into different bins. The absolute uncertainties in
stellar masses are expected to be large, but the relative uncertainties should
be small, so ignore them.

##  3 Alternative data treatment and tests

###  3.1 Data combination in the image plane

In addition to our $uv$-combined images, we also tried combining archival ALMA
data-product images (each one being produced individually using tclean in MFS
mode) in real space. This provides a consistency check for how well we are
able to combine in $uv$ space the data spanning a wide frequency range and
large array configuration range. For this test we focus on the deep central
mosaic (following the ASPECS footprint), where the data sets included are
listed in Table 1.

The first step in this process is to convolve each image to the same
resolution. A common beamsize was selected as the largest beam across all of
the data sets, and for each image the convolution kernel required to produce
this common beam was found numerically using the create_matching_kernel
function from the python photutils module (Bradley et al., 2022). Each image
was then convolved with this kernel. We expect that convolution to a lower
resolution will effectively decrease the sensitivity to point sources. The
largest beam across all of our data comes from the ASPECS pilot programme
(2013.1.00718.S, see Aravena et al. 2016), which has a beam of about
$1.5\,{\rm arcsec}\,{\times}\,0.8\,{\rm arcsec}$. The most sensitive map
(across a significantly larger area) comes from the full ASPECS survey
(2016.1.00324.L, see González-López et al. 2020), where the beam is about
$1.2\,{\rm arcsec}\,{\times}\,0.8\,{\rm arcsec}$, or the second-largest beam
across all of our images. We found that degrading the resolution of the full
ASPECS map to match the ASPECS pilot map ultimately produces a map with a
higher rms (after combining all the data) compared to removing the ASPECS
pilot map from the combination and convolving each image to match that of the
full ASPECS map; our final combined image therefore does not contain data from
the ASPECS pilot, although this amounts to only a small loss of data (see
Table 1).

We next created a pixel grid for the image. The size of the pixels in the grid
was chosen as the largest pixel size within the set of 61 images (in this case
0.2 arcsec), and the grid was set to span 5 arcmin, centred at 03h32m39.0s
$-$27∘47′29.1′′. Next, each of the 61 convolved images (minus those from the
ASPECS pilot) were reprojected onto this grid using the python function
reproject_interp, part of the reproject module (Robitaille et al., 2020).

The next step was to create a noise map for each reprojected and convolved
image. To do this, we created primary-beam-uncorrected images (simply the
primary beam image multiplied by the primary beam), and for each image we
created a mask using the same method described above. We then calculated the
rms of the primary-beam-uncorrected map multiplied by the mask, and divided
this by the primary beam.

These images should now be aligned to the same pixel grid, have the same
resolution, and have associated noise maps. The final step is to combine them,
weighting each pixel by its inverse variance. In order to remove boundary
artefacts associated with combining images with very different noise levels,
an apodizing mask was applied to the edges of each image in the combination.
The apodizing function chosen was a Gaussian with a standard deviation of the
beamsize, and the apodization was applied out to 3 times the beamsize from the
edge of each map. We investigated various options for apodization (e.g. a
cosine function or different apodization lengths), and found that this simple
choice removed most of the obvious edge effects.

Lastly, a new noise map was calculated for the combined image using the same
algorithm used to calculate a noise map for the $uv$-combined image. For
reference, the image-combined map can also be
downloaded.555https://doi.org/10.5683/SP3/YWBVWH

Figure 1: Left: Signal map of the deepest region of the HUDF after combining
all of the available archival ALMA Band-6 data given in Table 1. This region
is defined by the contour where the primary beam reaches 0.2, covering a total
of 4.2 arcmin2. Middle: Corresponding noise map. Right: S/N map, from dividing
the signal map by the noise map.

###  3.2 Comparison between uv-plane combination and image-plane combination

The optimal way to combine interferometric images is to add the observations
in the $uv$ plane, then image the entire data set. However, the different
frequencies and $uv$ sampling of the various observations may lead to unwanted
behaviour. Moreover, adding the individual pointings is much more
straightforward in the image plane than in $uv$ space. Therefore we would like
to compare the properties of the two images ($uv$-space combination and real-
space combination) to check that they are consistent with one another, and to
see which one performs better.

We focus on the region defined by our deep central map with a 250 k$\lambda$
taper, which was made going out to 0.2 times the primary beam (see Fig. 1). We
extracted all sources with S/N$\,{>}\,4$ within this region from both maps
using the same source-extraction procedure as described in Section 4; this
threshold was simply chosen so that enough overlapping sources could be
extracted for comparison. There are 36 peaks with this significance in the
$uv$-space combined map and 29 peaks in the image-space combined map. In Fig.
2 we show all of the sources found with S/N$\,{>}\,4$, from which it can be
seen that all bright and obvious sources are in agreement. However, we do see
more sources detected in the $uv$ combination compared to the image
combination, especially around the edges of our selected region, where the
$uv$ map does much better.

To quantify the difference between the two maps, in Fig. 3 we plot the ratio
of the peak S/N from the $uv$ combination to the peak S/N from the image
combination as a function of 1-mm flux density in the $uv$-space combined map.
For the galaxies with S/N$\,{>}\,4$ in the $uv$-combined map but not the
image-combined map, we simply extract the value of the image-combined map at
the location of the $uv$ detections in order to include them on this plot. We
find a relatively uniform improvement in S/N from combining the $uv$ data as a
function of source brightness. The mean S/N of all matching sources is higher
by about 1.2. Of course there is a balance between the increased complexity
and computing resources required to combine the data in the $uv$ plane versus
combining in the image plane, but our results indicate that the improvement
from combining the data in the $uv$-plane is worthwhile.

Figure 2: Results from combining all of the data in Table 1 in the $uv$ plane
before imaging, and imaging the same data individually before combining.
Signal maps are shown in the two top left panels, noise maps are shown in the
two top right panels and S/N maps are shown in the bottom two panels. In each
case, the $uv$ combination is on the left and the real-space image combination
is on the right. Sources with S/N$\,{>}\,$4 are marked with green squares in
the S/N maps. Figure 3: Ratio of the peak S/N from the $uv$ combination to the
peak S/N from the image combination, as a function of 1230-$\mu$m flux density
in the image combination. The peak S/N values from the $uv$ combination are
generally higher than for the image combination (the mean ratio is 1.2).

##  4 Results for sources

###  4.1 Cataloguing 1-mm sources

The central deep map shown in Fig. 1 contains all existing ALMA Band-6 data in
this region and thus should be the best map currently achievable; for
reference, the deepest part goes down to 4.6 $\mu$Jy beam-1, which is 50 per
cent deeper than any previous map. Over an area of 1.5 arcmin2 our new map has
an rms below 9.0 $\mu$Jy beam-1, which is 5 per cent better than the deepest
previously available map of this size (from ASPECS), and it has a noise level
below 35 $\mu$Jy beam-1 over the full $4.2\,{\rm arcmin}^{2}$. It is therefore
of interest to see if any new galaxies are detected in this new map.

We ran the simple peak-finding algorithm find_peaks, available in the
photutils python module, on both the S/N map contained within the region of
interest, and on the negative of the same S/N map. One approach to setting the
detection threshold is to use the most significant negative peak to set the
level above which we might expect all positive peaks to be real galaxies; for
reference, the most negative peak was found to have a value of $-4.1\,\sigma$,
and there are 35 positive peaks brighter than $+4.1\,\sigma$. However, we can
also lower the threshold to include more real sources at the expense of being
less confident about the reality of each one.

As an alternative means of setting the threshold, we looked at how the ratio
of the total number of positive to negative peaks greater than a given S/N
thresholds varies. We can define the ‘purity’ to be $1\,{-}\,{\cal N}_{\rm
neg}/{\cal N}_{\rm pos}$, such that a value of 0 means there are as many
negative as positive peaks and 1 is reached when there are no more significant
negative peaks. In Fig. 4 we plot this purity as a function of S/N threshold.
The choice of a threshold is of course a trade off (a smaller threshold will
result in more false positives), and to include more candidate sources (with
the understanding that not all may be real) we choose a purity of 0.7, which
corresponds to ${\rm S/N}=3.6$. There are 45 sources above this threshold. The
catalogue from ASPECS used a ‘fidelity’ threshold (the differential ratio of
galaxies in S/N bins) that resulted in a least significant source with ${\rm
S/N}\,{=}\,3.3$, while Dunlop et al. (2017) used a S/N threshold of 3.5 but
removed sources with no HST counterpart. Our choice of threshold is therefore
similar to (although slightly more conservative than) what was used in
previous studies, and we can also use our JADES catalogue to investigate
possible false positives.

Figure 4: Purity, i.e. 1 minus the the total number of positive to negative
peaks greater than a given S/N map, plotted as a function of S/N. A purity of
0.7 corresponds to a S/N threshold of 3.6, which we use to extract sources
from our deep map (Fig. 1).

The positions of sources extracted from this search are shown in Fig. 6, with
positions and 1-mm flux densities given in Table 2. Appendix B contains 5
arcsec$\,{\times}\,$5 arcsec cutouts overlaid on JWST F356W images. For
single-dish surveys, DSFGs are almost always unresolved, but this is not the
case for ALMA data. Hence we need to decide how to quote brightness values
when some DSFGs are resolved. For the flux densities, we follow a procedure
similar to the ASAGAO survey; for each source we fit a 2-dimensional Gaussian
profile, fixing the position to the peak pixel and the amplitude to the value
of the peak pixel, but allowing the size, ellipticity and position angle to
vary. We then calculate the number of beams contained within each source as

${\cal N}_{\mathrm{b}}=\frac{ab}{\theta_{\mathrm{maj}}\theta_{\mathrm{min}}},$
(1)

where $a$ and $b$ are the best-fit major and minor FWHM values, and
$\theta_{\mathrm{maj}}$ and $\theta_{\mathrm{min}}$ are the synthesized beam
major and minor FWHM, respectively. If the number of beams is greater than 1,
we calculate the integrated flux density at 243 GHz, or $1230\,\mu$m, as

$S_{1230}=S_{\mathrm{peak}}{\cal N}_{\mathrm{b}},$ (2)

otherwise the integrated flux density is simply the peak pixel value. We also
flag fits where the minor/major axis ratio is less than 0.5 as bad fits, and
use peak pixel values for these sources. Uncertainties are taken from the
noise map, and uncertainties from the fits are propagated to the integrated
flux densities. All of our sources and 1-mm flux densities are listed in Table
2.

As a check, we compare the 1-mm flux densities extracted here to the flux
densities given in the four previously published surveys. For simplicity, in
this comparison we simply match published sources to our new catalogue using a
search radius of 1 arcsec. The mean frequency of the map from Dunlop et al.
(2017) is 221 GHz (see Table 1), so we correct their flux densities to the
mean frequency of our map (243 GHz) assuming a modified blackbody SED with
spectral index $\beta\,{=}\,1.5$, a dust temperature of 30 K, and a redshift
of 1.5; the correction factor is 1.3. Similarly, the mean frequency of the map
from GOODS-ALMA (Gómez-Guijarro et al., 2022) is 265 GHz, so we follow the
same procedure and apply a correction factor of 0.8. The mean frequency of the
remaining maps are effectively identical to ours, so we do not apply any
further corrections.

The flux densities of matched peaks are shown in Fig. 5, and the cross-matched
IDs are given in Table 3. We find generally good agreement with all of the
previously-published flux densities after applying the above corrections.

In Fig. 6 we show our detected sources compared to those found by Dunlop et
al. (2017), the ASAGAO survey (Hatsukade et al., 2018) and the ASPECS survey
(González-López et al., 2020), as well as a few sources from the wider GOODS-
ALMA survey (Gómez-Guijarro et al., 2022). It appears that most of the
detections in our combined map coincide with published sources, but there are
13 new DSFGs. There are also a few published detections that do not appear
here; typically, these are low-significance sources that could have been
positive noise excursions.

In Appendix A we perform the same source extraction procedure on the larger
but shallower map covering the ASAGAO footprint (excluding the central deep
region, where we have already detected all of the sources in this shallow
map), now fixing the detection threshold to $4.5\,\sigma$ (as was done to make
the ASAGAO catalogue). We find that the purity at this threshold is 0.9 (which
is therefore fairly conservative) and we find a total of 27 additional
galaxies, nine of which are new. We use a similar algorithm to measure flux
densities, and find similar flux densities compared to those published by
ASAGAO (Hatsukade et al., 2018) and GOODS-ALMA (Gómez-Guijarro et al., 2022).
As a test, we also extract peaks from the central region of the shallow map
and measured their flux densities, and find good agreement with the flux
densities measured in our deep map.

Figure 5: Comparison of the flux densities measured in this work with flux densities from Dunlop et al. (2017), ASPECS (González-López et al., 2020), ASAGAO (Hatsukade et al., 2018) and GOODS-ALMA (Gómez-Guijarro et al., 2022). Figure 6: Signal-to-noise ratio map of the deep central region, covering $4.2\,{\rm arcmin}^{2}$, with peaks ${>}\,3.6\,\sigma$ indicated as green boxes. Galaxies found by Dunlop et al. (2017) are shown as yellow circles, galaxies from the ASPECS survey (González-López et al., 2020) are orange circles and galaxies from the ASAGAO survey (Hatsukade et al., 2018) are red circles. A few of the brightest galaxies are also detected in the GOODS-ALMA survey (Gómez-Guijarro et al., 2022) and shown as brown circles. The gold contour shows the footprint of the deepest HST data from the CANDELS/HUDF09 survey, while the entire region is covered by the JADES survey. Galaxies with counterparts in both the CANDELS catalogue and the JADES catalogue are indicated with a blue cross and galaxies with only a JADES counterpart are indicated with a red cross. There are no galaxies with a CANDELS/HUDF09 counterpart but no JADES counterpart. Numbers refer to the labels in Tables 2 and 3. Table 2: Positions and flux densities ($S_{1230}$) for the 45 sources found in the ALMA 1-mm combined image at ${\rm S/N}\,{>}\,3.6$. Here flux densities are measured by fitting Gaussian profiles to the sources in order to estimate the number of beams per source, then scaling the peak pixel value by the number of beams; for sources consistent with 1 beam or less, we set the flux density to be the peak pixel value. For galaxies with matches in the JADES catalogue and sufficient photometry to fit SEDs (S/N$\,{>}\,$5 in the F356W band) we provide stellar masses and photometric redshifts from McLeod et al. (in prep.). Galaxies with spectroscopic redshifts from Rieke & the JADES Collaboration (2023) are indicated by asterisks, which we provide instead of photometric redshifts. The median photometric redshift uncertainty of our ALMA sample is $\Delta z/(1\,{+}\,z)\,{\simeq}\,0.03$, while the relative uncertainties in stellar masses are small. Corresponding identifications in other published catalogues are given in Table 3. ALMA ID 45 is at the position of a low-$z$ and a high-$z$ galaxy. We assume that this ALMA source corresponds to the high-$z$ galaxy, the photometric redshift and stellar mass are derived using only long-wavelength JWST photometry where the high-$z$ galaxy is detected. Name | RA Dec [J2000] |  $S_{1230}$ [$\mu$Jy] | $z$ | $\log\left(M_{\ast}/\mathrm{M}_{\odot}\right)$
---|---|---|---|---
ALMA-HUDF-1 | 3:32:43.53 $-$27:46:39.2 | 850$\pm$23 | 2.850∗ | 11.0
ALMA-HUDF-2 | 3:32:38.54 $-$27:46:34.6 | 749$\pm$10 | 2.543∗ |  09.8
ALMA-HUDF-3 | 3:32:36.96 $-$27:47:27.2 | 507$\pm$12 | 2.470∗ | 10.7
ALMA-HUDF-4 | 3:32:39.75 $-$27:46:11.6 | 488$\pm$21 | 1.546∗ | 11.1
ALMA-HUDF-5 | 3:32:34.44 $-$27:46:59.8 | 477$\pm$17 | 1.413∗ | 10.7
ALMA-HUDF-6 | 3:32:40.08 $-$27:47:55.6 | 381$\pm$25 | 1.998∗ | 10.7
ALMA-HUDF-7 | 3:32:41.01 $-$27:46:31.6 | 342$\pm$13 | 2.400∗ | 10.4
ALMA-HUDF-8 | 3:32:43.33 $-$27:46:47.0 | 258$\pm$12 | 2.410∗ | 10.6
ALMA-HUDF-9 | 3:32:35.08 $-$27:46:47.8 | 231$\pm$12 | 2.497∗ | 10.8
ALMA-HUDF-10 | 3:32:35.56 $-$27:47:04.2 | 189$\pm$12 | 3.450∗ |  09.9
ALMA-HUDF-11 | 3:32:39.88 $-$27:47:15.2 | 181$\pm$19 | 1.095∗ | 10.5
ALMA-HUDF-12 | 3:32:38.03 $-$27:46:26.6 | 168$\pm$90 | 1.159∗ |  09.6
ALMA-HUDF-13 | 3:32:37.35 $-$27:46:45.8 | 161$\pm$20 | 1.846∗ |  09.9
ALMA-HUDF-14 | 3:32:35.51 $-$27:46:26.8 | 154$\pm$36 | 1.070∗ | 10.1
ALMA-HUDF-15 | 3:32:42.99 $-$27:46:50.2 | 134$\pm$90 | 1.036∗ | 10.7
ALMA-HUDF-16 | 3:32:36.48 $-$27:46:31.8 | 133$\pm$14 | 1.120∗ |  09.5
ALMA-HUDF-17 | 3:32:38.80 $-$27:47:14.8 | 130$\pm$16 | 1.850∗ | 10.5
ALMA-HUDF-18 | 3:32:42.37 $-$27:47:07.8 | 126$\pm$13 | 1.317∗ | 11.0
ALMA-HUDF-19 | 3:32:36.19 $-$27:46:28.0 | 119$\pm$14 | 2.307∗ | 10.6
ALMA-HUDF-20 | 3:32:38.75 $-$27:48:10.4 | 114$\pm$19 | 2.940∗ | 10.1
ALMA-HUDF-21 | 3:32:41.69 $-$27:46:55.6 |  098$\pm$80 | 1.776∗ | 10.3
ALMA-HUDF-22 | 3:32:34.86 $-$27:46:40.8 |  096$\pm$19 | 1.098∗ | 10.2
ALMA-HUDF-23 | 3:32:41.04 $-$27:47:48.0 |  089$\pm$22 | … |  0…
ALMA-HUDF-24 | 3:32:35.78 $-$27:46:27.6 |  086$\pm$16 | 1.094∗ | 10.4
ALMA-HUDF-25 | 3:32:41.24 $-$27:46:16.6 |  081$\pm$22 | … |  0…
ALMA-HUDF-26 | 3:32:34.71 $-$27:46:45.0 |  078$\pm$22 | 1.450∗ | 10.1
ALMA-HUDF-27 | 3:32:35.99 $-$27:47:25.8 |  077$\pm$15 | 2.580∗ | 10.3
ALMA-HUDF-28 | 3:32:34.82 $-$27:46:31.2 |  076$\pm$21 | … |  0…
ALMA-HUDF-29 | 3:32:38.50 $-$27:47:02.6 |  066$\pm$15 | 0.954∗ | 10.9
ALMA-HUDF-30 | 3:32:37.61 $-$27:47:44.0 |  066$\pm$10 | 1.240∗ |  09.5
ALMA-HUDF-31 | 3:32:41.45 $-$27:47:29.2 |  058$\pm$14 | 0.621∗ |  09.0
ALMA-HUDF-32 | 3:32:37.73 $-$27:47:07.2 |  057$\pm$14 | 0.668∗ |  09.8
ALMA-HUDF-33 | 3:32:38.59 $-$27:47:30.4 |  054$\pm$14 | 2.642∗ |  09.6
ALMA-HUDF-34 | 3:32:40.23 $-$27:47:38.2 |  050$\pm$11 | … |  0…
ALMA-HUDF-35 | 3:32:35.75 $-$27:46:39.4 |  046$\pm$14 | 2.070∗ |  09.8
ALMA-HUDF-36 | 3:32:35.77 $-$27:46:55.4 |  042$\pm$12 | 1.721∗ |  09.8
ALMA-HUDF-37 | 3:32:38.56 $-$27:47:05.4 |  041$\pm$90 | … |  0…
ALMA-HUDF-38 | 3:32:41.83 $-$27:46:56.8 |  038$\pm$80 | … |  0…
ALMA-HUDF-39 | 3:32:39.78 $-$27:46:29.4 |  038$\pm$11 | … |  0…
ALMA-HUDF-40 | 3:32:37.08 $-$27:46:17.4 |  038$\pm$70 | 2.226∗ |  09.0
ALMA-HUDF-41 | 3:32:38.69 $-$27:46:30.8 |  034$\pm$90 | … |  0…
ALMA-HUDF-42 | 3:32:43.62 $-$27:46:59.0 |  033$\pm$80 | 1.569∗ |  09.7
ALMA-HUDF-43 | 3:32:36.37 $-$27:46:50.2 |  032$\pm$80 | … |  0…
ALMA-HUDF-44 | 3:32:41.90 $-$27:46:58.6 |  032$\pm$80 | … |  0…
ALMA-HUDF-45 | 3:32:42.35 $-$27:46:57.4 |  030$\pm$70 | 2.420∗ |  09.8
Table 3: Cross-matched identifications between our catalogue of sources and previous studies. Near-infrared matches are from JADES (Rieke & the JADES Collaboration, 2023) and CANDELS (Guo et al., 2013). Matches with other ALMA Band-6 surveys are from Dunlop et al. (2017), ASAGAO (Hatsukade et al., 2018), GOODS-ALMA (Gómez-Guijarro et al., 2022), and ASPECS (González-López et al., 2020). ALMA ID 44 is not detected in our F356W-selected catalogue, but has a counterpart in the JADES catalogue, which we provide below. The JADES ID corresponding to ALMA ID 45 is a blend of a low-$z$ and a high-$z$ galaxy, and we assume that our ALMA source corresponds to the galaxy at high-$z$. Name | JADES ID | CANDELS ID | Dunlop et al. ID | ASAGAO ID | GOODS-ALMA ID | ASPECS ID
---|---|---|---|---|---|---
ALMA-HUDF-1 | 208820 | J033243.52$-$274639.0 | UDF2 | 4 | A2GS9 | C06
ALMA-HUDF-2 | 209117 | J033238.54$-$274634.0 | UDF3 | 5 | A2GS25 | C01
ALMA-HUDF-3 | 204232 | J033236.96$-$274727.2 | UDF5 | 12 | A2GS41 | C02
ALMA-HUDF-4 | 211273 | J033239.73$-$274611.2 | UDF8 | 16 | … | C05
ALMA-HUDF-5 | 207012 | J033234.43$-$274659.5 | UDF6 | 13 | … | C03
ALMA-HUDF-6 | 202563 | J033240.05$-$274755.4 | UDF11 | 15 | … | C10
ALMA-HUDF-7 | 209357 | J033241.02$-$274631.4 | UDF4 | 10 | … | C04
ALMA-HUDF-8 | 208030 | J033243.32$-$274646.7 | UDF7 | 14 | … | C11
ALMA-HUDF-9 | 208000 | J033235.07$-$274647.5 | UDF13 | 23 | … | C07
ALMA-HUDF-10 | 206834 | J033235.55$-$274703.8 | … | … | … | C09
ALMA-HUDF-11 | 205449 | J033239.88$-$274715.0 | … | … | … | C16
ALMA-HUDF-12 | 209777 | J033238.02$-$274626.2 | … | … | … | C08
ALMA-HUDF-13 | 208134 | J033237.35$-$274645.4 | … | … | … | C18
ALMA-HUDF-14 | 209492 | J033235.48$-$274626.6 | … | … | … | C23
ALMA-HUDF-15 | 207739 | J033242.98$-$274649.9 | … | … | … | C13
ALMA-HUDF-16 | 209285 | … | … | … | … | C12
ALMA-HUDF-17 | 205379 | J033238.79$-$274714.7 | … | … | … | C17
ALMA-HUDF-18 | 206183 | J033242.37$-$274707.6 | UDF16 | … | … | C15
ALMA-HUDF-19 | 209617 | J033236.17$-$274627.6 | … | … | … | C19
ALMA-HUDF-20 | 201501 | … | … | … | … | C24
ALMA-HUDF-21 | 207227 | J033241.68$-$274655.4 | … | … | … | C14a
ALMA-HUDF-22 | 208277 | J033234.85$-$274640.4 | … | … | … | C25
ALMA-HUDF-23 | … | … | … | … | … | …
ALMA-HUDF-24 | 209480 | J033235.77$-$274627.4 | … | … | … | C20
ALMA-HUDF-25 | … | … | … | … | … | …
ALMA-HUDF-26 | 208267 | … | … | … | … | C26
ALMA-HUDF-27 | 204579 | J033235.98$-$274725.6 | … | … | … | C21
ALMA-HUDF-28 | 129574 | … | … | … | … | …
ALMA-HUDF-29 | 206703 | J033238.48$-$274702.4 | … | … | … | C33
ALMA-HUDF-30 | 203384 | J033237.61$-$274744.0 | … | … | … | C22
ALMA-HUDF-31 | 204483 | J033241.45$-$274729.3 | … | … | … | …
ALMA-HUDF-32 | 206205 | J033237.73$-$274706.9 | … | … | … | C32
ALMA-HUDF-33 | 204449 | J033238.55$-$274730.2 | … | … | … | …
ALMA-HUDF-34 | … | … | … | … | … | C27
ALMA-HUDF-35 | 208812 | J033235.73$-$274639.0 | … | … | … | …
ALMA-HUDF-36 | 124908 | J033235.76$-$274655.0 | UDF15 | … | … | …
ALMA-HUDF-37 | … | … | … | … | … | …
ALMA-HUDF-38 | … | … | … | … | … | C14b
ALMA-HUDF-39 | 265959 | J033239.76$-$274629.3 | … | … | … | …
ALMA-HUDF-40 | 210730 | J033237.07$-$274617.1 | … | … | … | C31
ALMA-HUDF-41 | … | … | … | … | … | …
ALMA-HUDF-42 | 207079 | J033243.61$-$274658.7 | … | … | … | …
ALMA-HUDF-43 | … | … | … | … | … | …
ALMA-HUDF-44 | 267661 | J033241.92$-$274658.2 | … | … | … | …
ALMA-HUDF-45 | 207277 | J033242.35$-$274657.0 | … | … | … | …

###  4.2 Matching mm-selected sources to near-IR-selected sources

The expected uncertainty in our ALMA positions is
$\delta\mathrm{RA}=\delta\mathrm{Dec}=0.6\,{\times}\,{\rm
FWHM}\,{\div}\,{\rm(S/N)}$ (Ivison et al., 2007), so we expect the $1\,\sigma$
radial uncertainty to be $\delta r\,{\approx}\,$0.76 arcsec/(S/N) (using the
geometric mean of the elliptical ALMA beam). Since the probability density of
finding a source a distance $r$ from its true position is proportional to
$r\mathrm{e}^{-r^{2}/2\delta r^{2}}$, one must go out to a distance of
2.5$\delta r$ in order to find a correct match with 95 per cent certainty. For
a given ALMA detection, we thus search the JADES catalogue out to a distance
of 1.9 arcsec/(S/N), where S/N here is the S/N of each source found above our
detection threshold. For the most significant ALMA sources this search radius
is unphysically small (much less than 1 pixel), so we apply a minimal search
radius of 1.5 pixels.

We preform the same counterpart search with the CANDELS catalogue, including
the deep infrared data from the HUDF09 survey (Bouwens et al., 2011). It
should be noted that Dunlop et al. (2017) found a systematic offset between
the CANDELS and ALMA astrometry, which was resolved by applying an offset of
about 0.25 arcsec south to the CANDELS positions; the same offset was applied
here. In Table 3 we provide the JADES and CANDELS IDs for these matches, and
in Appendix B we show the positions of matched sources in our ALMA cutouts
overlaid over JWST F356W imaging.

We note that one ALMA source (ID 44) did not have a counterpart in our
F356W-selected catalogue yet had a counterpart in the JADES catalogue, so we
include this ID here. At the position of our ALMA source ID 45 we found that
there were two blended galaxies in the longer-wavelength JWST images, yet at
shorter JWST wavelengths one of these galacxies was undetected. In the JADES
catalogue there was only one ID for these two sources (ID 207277) that had a
spectroscopic redshift of 0.332. This spectroscopic redshift likely
corresponds to the galaxy brighter at short wavelengths, while our ALMA source
is probably the second galaxy that is only detected at longer wavelengths. We
therefore use the longer wavelength photometry (where this galaxy is detected)
to fit an SED, which results in a photometric redshift of 2.42. Throughout the
rest of the paper we report results derived from this fit.

Sources found to have counterparts within their given search radii in both the
CANDELS catalogue and the JADES catalogue are marked in Fig. 6 with blue
crosses. There are 35 ALMA-detected galaxies with counterparts in both
catalogues, seven of which are ALMA sources that have not been previously
reported. There are an additional four ALMA-detected galaxies with a
counterpart only in the JADES catalogue, one of which is an ALMA source that
has not been previously reported. We do not find any CANDELS counterparts with
no corresponding JADES counterpart, which highlights one of the benefits of
the deeper catalogue of the HUDF made possible with JWST (although the
improvement is not dramatic, since CANDELS gets close to identifying all of
our mm sources).

This leaves seven ALMA-detected sources with no optical or near-infrared
counterpart; two of these have previously-published ALMA identifications, the
remaining five do not. It is worth pointing out that two of the sources with
no counterparts are close to the edge of the map where the noise is the
largest (IDs 23 and 25), but the other five are found near the centre of the
map where the data are much deeper (IDs 34, 37, 38, 41 and 43).

The number of beams in our ALMA map can be estimated as the area of the map
(4.2 arcmin2) divided by the beam area ($2\pi\theta_{\rm maj}\theta_{\rm
min}/8\ln 2$), which results in about 8500. From Gaussian statistics, we would
expect around 1.5 beams to randomly have a significance greater than
$3.6\,\sigma$; however, there are really twice as many statistically
independent noise samples due to correlations with the beam (e.g., Condon,
1997; Condon et al., 1998; Dunlop et al., 2017), so really we would expect
three false positive detections. The measured number of negative peaks more
significant than 3.6 is 14, which is larger than the expectation from a
Gaussian distribution – the negative peaks might not be perfectly Gaussian
distributed, but this is not surprising given the non-uniform $uv$ coverage.
On the other hand, the number of ALMA-detected sources with no optical or
near-infrared counterpart (seven) is more comparable to the number of negative
peaks, and so we cannot rule out the possibility that some of them are false
positives.

In Fig. 7 we show the distributions of redshifts, magnitudes, and colours from
our JADES catalogue. For the redshifts and colours, we filter out three
galaxies which have a S/N$\,{<}\,$5 in the F356W band as the photometry is not
reliable. For the remaining galaxies, the median redshift uncertainty is
$\Delta z/(1\,{+}\,z)\,{\simeq}\,0.03$, so we expect them to be accurate. We
display the distributions of all 38/45 ALMA galaxies with JWST detections and
also highlight the distributions of the eight ALMA galaxies that were not
found in previous surveys. The first (top left) panel in Fig. 7 shows the
distribution of redshifts (spectroscopic where available, otherwise
photometric). We see that the redshift distribution is flat around
$z\,{\simeq}\,2$, and our new galaxies appear to follow this trend. Next we
show the distributions of the magnitudes of the sources in the NIRCam images
(F277W, F335M, F356W, F410M and F444W). We see that most ALMA sources have
magnitudes ranging from 26 to 32, with the three faintest sources extending
out to 40 (these are the three sources that we do not fit SEDs to). Lastly we
show distributions for three NIRCam colours; most ALMA galaxies are fairly red
(i.e. brighter at longer wavelengths), and our new, fainter ALMA galaxies tend
to follow the same distribution.

Figure 7: Distribution of ALMA source properties matched to our JADES
catalogue (Rieke & the JADES Collaboration 2023, McLeod et al. in prep.). The
first (top left) panel shows the distribution of redshifts (spectroscopic when
available, otherwise photometric). The next five panels show distributions of
NIRCam magnitudes for the ALMA sources in the JADES catalogue. The final three
panels show distributions of three selected NIRCam colours. The distributions
from all 38/45 ALMA galaxies with JWST detections are shown in blue, while the
distributions from the subset of eight new ALMA galaxies with JWST detections
are shown in red. In the redshift panel (top left) there are only 35 ALMA
galaxies with available redshifts, and five new ALMA galaxies with redshifts.
In the bottom row we also only plot those galaxies with redshifts because the
colours are quite uncertain.

###  4.3 Stellar mass-redshift distribution

As described in Section 2.4.1, the photometric redshifts and stellar masses
for most JADES galaxies have been estimated using the extensive
multiwavelength imaging available for the HUDF. It is therefore of interest to
see how the photometric redshifts and stellar masses of our mm-selected
galaxies compare to the typical galaxies found in this field. For our list of
sources with a JADES counterpart, we filter out all sources with S/N$\,{<}\,$5
in the F356W band, since these sources will not have reliable SED fits. For
the remaining sources, we check to see if there is a spectroscopic redshift,
and otherwise use the photometric redshift.

In Fig. 8 we plot stellar mass versus redshift for all JADES galaxies (McLeod
et al. in prep.), and highlight our ALMA Band-6-selected sources with good SED
fits in red. We find that most of our ALMA-selected sources are high-stellar-
mass galaxies between $z\,{=}\,1$ and 3, similar to what was found in earlier
ALMA surveys of the HUDF (e.g., Dunlop et al., 2017; McLure et al., 2018;
Aravena et al., 2020). In particular, there are 53 galaxies with
$M_{\ast}\,{>}\,10^{10}\,$M⊙, 21 of which have been selected in our ALMA
image. Also of note is that at 1.23 mm we are sensitive primarily to
$z\,{>}\,1$ objects, since in our sample of 35 galaxies with spectroscopic or
photometric redshifts, only three are at $z\,{<}\,1$. In a similar vein, it
may also be worth noting that all of the $\log(M_{\ast}/{\rm
M}_{\odot})\,{>}\,10.3$ galaxies at $2\,{<}\,z\,{<}\,3$ are detected in our
ALMA image.

In order to investigate the difference between galaxies with
$M_{\ast}\,{>}\,10^{10}\,$M⊙ that we have detected with ALMA versus those that
we have not detected with ALMA, in Fig. 9 we show the distributions of three
select NIRCam magnitudes and three NIRCam colours for the subsamples of 21
ALMA-detected galaxies and 32 ALMA-undetected galaxies. We find that there is
no discernible difference in magnitudes between the two samples; however, mm-
bright ALMA sources tend to have redder NIRCam colours than the galaxies with
similarly high stellar masses that we have not detected with ALMA.

Figure 8: Stellar mass versus redshift for JADES galaxies (McLeod et al. in
prep.), with JADES galaxies detected in our $S_{1230}$ image highlighted in
red (circles indicate photometric redshifts, stars indicate spectroscopic
redshifts). The apparent vertical features are a result of the grid used for
the photometric redshift fitting. Figure 9: Top row: Distributions of three
select NIRCam magnitudes for all ALMA-detected galaxies with JWST counterparts
and with stellar masses above 10${}^{10}\,$M⊙, compared to all JWST galaxies
with stellar masses above 10${}^{10}\,$M⊙ that are not detected by ALMA.
Bottom row: Same as top row but for three select NIRCam colours. While we do
not see a significant difference between the two samples in terms of
magnitude, the mm-bright galaxies detected by ALMA with high stellar masses
tend to have redder NIRCam colours compared to ALMA-undetected galaxies with
equally large stellar masses.

###  4.4 Stacking on near-IR-selected positions

One major benefit of the deep ancillary catalogues available in the HUDF is
the ability to stack on the positions of undetected galaxies in different
stellar mass and redshift bins, thus estimating the statistical properties of
fainter mm-emitting sources that are not individually detectable in our image.

To do this, we follow a procedure similar to Simstack (Viero et al., 2013b).
However, since we adapt this for a non-circular beam and non-uniform noise
distribution, it is worth describing what we do in a little detail. First, we
mask out all galaxies that we have detected in our ALMA map (see Table 2),
with a source mask set to be $3\times{\rm FWHM}$ in diameter. We also restrict
this stacking analysis to the deep central map shown in Fig. 6. Then we
subtract the weighted mean of the image – this is crucial, since the ‘stack’
is really the covariance between a map and catalogue (see section 3 of Marsden
et al., 2009) and will give a biased result unless the map has zero mean.

Next, we define a grid of four stellar mass bins, logarithmically spaced
between $\log\left(M_{\ast}/\mathrm{M}_{\odot}\right)\,{=}\,7.4$ and
$\log\left(M_{\ast}/\mathrm{M}_{\odot}\right)\,{=}\,11.4$, with a width of
$\Delta\log\left(M_{\ast}/\mathrm{M}_{\odot}\right)\,{=}\,1$. We also define a
grid of four redshift bins between 0 and 8, with a width of 2. For the
stacking catalogue, we again turn to JADES, using stellar masses and redshifts
given by McLeod et al. (in prep.); here redshifts are spectroscopic if
available, otherwise photometric. These bins are chosen because they are
expected to contain the galaxies comprising the vast majority of the CIB (see
Section 5.3 for more details).

For each redshift and stellar mass bin we produce a ‘hits’ map, which is
simply a copy of our ALMA map with pixel values set to the number of JADES
galaxies within the bin contained within each pixel. We convolve each hits map
with the ALMA beam, thereby producing a model image of the sky where the only
free parameter of each map is its amplitude, which in this case can be
interpreted as the best-fit flux density of the galaxies within the given
stellar mass and redshift bin. As was done with the data, we also subtract the
weighted mean of the convolved maps.

As described in Viero et al. (2013b), in general there are correlations
between redshift bins simply due to the presence of large-scale structure, and
to take these into account one would have to fit for all of the amplitudes in
the defined bins simultaneously. In our case, redshift bins have a width of 2,
so these correlations are expected to be negligible. We therefore fit for the
amplitude in each bin independently, which reduces to solving a simple
weighted linear regression between the ALMA map and maps of ‘hits’ in the
JADES catalogue for each bin. The solution is

$\hat{S}_{\alpha}=\frac{\sum_{j}w_{j}D_{j}N^{j}_{\alpha}}{\sum_{j}w_{j}(N^{j}_{\alpha})^{2}},$
(3)

where $j$ labels each map pixel, $\alpha$ denotes the redshift/stellar mass
bin, $N_{\alpha}$ is the hits map of bin $\alpha$ convolved with the ALMA beam
with the weighted mean subtracted, $D$ is the data map with the weighted mean
subtracted and $w_{j}\,{=}\,1/\sigma_{j}^{2}$ are the inverse-variance
weights. The sum is performed over all the pixels in the map. Equation 3 is
effectively the covariance between the map and the catalogue, but with the
ALMA beam taken into account.

To estimate the uncertainties in $\hat{S}_{\alpha}$, for each redshift and
stellar mass bin we generate 1000 catalogues with the same number of sources
in each bin but with random positions, and evaluate Eq. 3 for each one. We
find the distribution of values to be well-described by a Gaussian, so we take
the standard deviation of the random catalogue flux densities to be the
$1\,\sigma$ error in $\hat{S}_{\alpha}$. In order to visually show our
results, we also evaluate Eq. 3 after shifting the hits map $N_{\alpha}$
relative to the data map $D$ within boxes of 10 arcsec; this is effectively
the cross-correlation between the two maps, where the central pixel is the
zero-lag cross-correlation (or the covariance), which is the value we are most
interested in.

The top panel of Figure 10 shows our results from the cross-correlation, while
in the bottom panel we list the central pixel values, which are the best-fit
flux densities of the galaxies within each bin. For discussing the CIB
contributions we want to know the pixel sum of the best-fit average surface
brightness for each bin (weighted by the inverse-variance), which we calculate
as

$\begin{split}\hat{I}_{\alpha}&=\frac{\hat{S}_{\alpha}}{2\pi\sigma_{\rm
maj}\sigma_{\rm min}}\frac{\sum_{j}w_{j}N_{\alpha}^{j}}{\sum_{j}w_{j}}\\\
&=\frac{\hat{S}_{\alpha}N_{\alpha,\mathrm{eff}}}{A}.\end{split}$ (4)

Here $\sigma_{\rm maj}$ and $\sigma_{\rm min}$ are the beam major and minor
axes (in standard deviation units), respectively, the quantity $A$ is the
solid angle of the map and $N_{\alpha,\mathrm{eff}}$ is the effective total
number of sources from catalogue $\alpha$ in the map, weighted by the noise,

$N_{\alpha,\mathrm{eff}}=\frac{\sum_{j}w_{j}N_{\alpha}^{j}}{\sum_{j}w_{j}}\frac{A}{2\pi\sigma_{\rm
maj}\sigma_{\rm min}}.$ (5)

This is just the average number of sources per beam (weighted by the noise),
times the number of beams in the map. In Fig. 10 we provide the values of
$N_{\alpha,\mathrm{eff}}$ along with $\hat{S}_{\alpha}$. Note that when we
calculate the weighted mean surface brightness of the sky from our model, the
weighted mean should not be subtracted from $N_{\alpha}^{j}$.

It is worth noting that a simpler stacking technique (just summing pixel
values at the positions of JWST galaxies; see Marsden et al. 2009) produces
similar results, merely with slightly less significance. As a completely
separate null test, we also stacked our ALMA map at the positions of 21
sources from the JADES catalogue flagged as stars that happen to fall within
our ALMA map. This stack is shown in Fig. 11 and we can see that it is
consistent with noise.

In Fig. 10 we see that many high stellar mass bins are blank; this is because
we have either detected all of the galaxies within these bins, or there were
no galaxies within these bins to begin with. We also see that there are
stacked peaks ${\gtrsim}\,3\,\sigma$ in all bins with stellar masses between
$10^{8.4}\,{\rm M}_{\odot}$ and $10^{10.4}\,{\rm M}_{\odot}$ and with
redshifts between 0 and 4. Across the lowest stellar mass the stacks tend to
be positive but with error bars overlapping 0, meaning that galaxies with
stellar masses ${<}\,10^{8.4}\,$M⊙ barely contribute to our ALMA map (a topic
that is explored further in Section 5.3).

In terms of stellar mass, the main conclusions of these stacking results are
that: (1) there are about 60 galaxies with stellar masses around
$10^{10}\,{\rm M}_{\odot}$ lying just below our $3.6\,\sigma$ ALMA threshold,
which have flux densities around 15 $\mu$Jy; (2) there are more than 300
galaxies with stellar masses around $10^{9}\,{\rm M}_{\odot}$ that can also be
statistically detected in the ALMA map, with individual flux densities of a
few $\mu$Jy; and (3) galaxies with stellar masses around $10^{8}\,{\rm
M}_{\odot}$ or below have flux densities that are too low to be detected in
the ALMA map, even statistically.

Figure 10: Top: Stacks at the positions of JADES galaxies in redshift and
stellar mass (from McLeod et al. in prep.) bins, after masking all of our
detected ALMA sources. The cutouts are 10 arcsec$\,{\times}\,10\,$arcsec. The
colour represents S/N, where the signal is the variance-weighted mean and the
noise is the uncertainty in the variance-weighted mean. Bottom: Central pixel
values and uncertainties of the top panel (i.e. our best estimate of the
average 1230-$\mu$m flux density of galaxies within each bin), with the number
of JADES galaxies contributing to each bin also given. For a definition of the
quantity $N_{\rm eff}$, see Eq. 5. Figure 11: Stack at the positions of 21
JADES sources flagged as stars that happen to lie within our ALMA 1-mm map; as
expected, we do not find a significant signal.

##  5 The cosmic infrared background in the HUDF

The absolute intensity of extragalactic light has been studied at all
wavelengths from $\gamma$-rays to the radio (Hill et al., 2018). Peaking at
around $200\,\mu$m, the CIB is the brightness of the extragalactic sky at
infrared wavelengths, averaged over the whole sky, after subtracting all
contributions originating from the Solar System and the Milky Way. This
absolute value tells us about the history of star formation and has been
measured using the FIRAS instrument onboard the COBE satellite (Fixsen et al.,
1998). The 1-mm background lies on the long-wavelength side of the CIB and we
can try to use our new map to estimate what fraction can be accounted for in
DSFGs.

###  5.1 Resolved source contribution to the cosmic infrared background

An important question is whether the intensity of the CIB can be recovered by
summing the contribution from known galaxies, or if there exists an additional
population of sources or a genuinely diffuse component of the CIB. This
question has been addressed by many previous studies at wavelengths around 1
mm (e.g. Penner et al., 2011; Viero et al., 2013b; Dunlop et al., 2017;
Hatsukade et al., 2018; González-López et al., 2020; Gómez-Guijarro et al.,
2022; Chen et al., 2023), with results typically in the tens of per cent
range, depending on the precise wavelength. Here we explore this question with
our new ALMA map at 1.23 mm and our new JADES catalogue. To start with, the
sum of our detected source flux densities (${>}\,3.6\,\sigma$) is
(7.33$\,{\pm}\,$0.10) mJy, and these sources are detected within an area of
$1.15\,{\times}\,10^{-3}\,{\rm deg}^{2}$. Therefore we have resolved a total
intensity of $(6.35\,{\pm}\,0.09)\,{\rm Jy}\,{\rm deg}^{-2}$ in individually-
detected DSFGs.

In addition to this, our stacking analysis (Fig. 10) demonstrates that our map
is also statistically sensitive to fainter galaxies. Multiplying
$\hat{S}_{\alpha}$ by $N_{\rm eff}$ in each bin shown in Fig. 10, and summing,
yields a flux density of $(2.6\,{\pm}\,0.5)$ mJy, or an intensity of
$(2.4\,{\pm}\,0.5)\,{\rm Jy}\,{\rm deg}^{-2}$ (where the area used in our
stack is $1.11\,{\times}\,10^{-3}\,{\rm deg}^{2}$, slightly less than the full
map due to our source mask); this stacking result is larger than has
previously been possible, due to the combination of a deeper ALMA map and
larger catalogue from JWST. Lastly, noting that there are no galaxies with
$S_{1230}\,{>}\,0.85$ mJy present in our map, we could also add a contribution
from brighter sources. The GOODS-ALMA survey (Gómez-Guijarro et al., 2022)
covered a much larger area with ALMA at 1 mm (about 72 arcmin2) and found 22
galaxies with $S_{1230}\,{>}\,0.85$ mJy, corresponding to an intensity of
$(1.31\,{\pm}\,0.02)\,{\rm Jy}\,{\rm deg}^{-2}$ (including the factor of 0.8
to convert their flux densities from 268 GHz to 243 GHz) and so we can add
this to our resolved CIB contribution as well.

In total, we directly or statistically detect a total CIB intensity of
$(8.7\,{\pm}\,0.5)\,{\rm Jy}\,{\rm deg}^{-2}$ in our map, and by adding in a
contribution from brighter sources, we estimate that the true value of the CIB
is $(10.0\,{\pm}\,0.5)\,{\rm Jy}\,{\rm deg}^{-2}$. We note that if we instead
perform our stacking analysis on the full map (without masking bright detected
sources) we obtain a CIB estimate of $(8.7\,{\pm}\,0.8)\,{\rm Jy}\,{\rm
deg}^{-2}$, consistent with our approach of detecting bright objects then
masking them to add the stacking result. Now we have to determine what
fraction of the background we have accounted for.

Figure 12: CIB intensity estimates at 1.23 mm (or 243 GHz). The stepped lines
show running sums of individual galaxy flux densities from this work (red),
Dunlop et al. (2017, brown) and ASPECS (González-López et al., 2020, teal).
The contribution to the CIB from $S_{1230}\,{>}\,0.85\,$mJy sources
($1.3\,{\rm Jy}\,{\rm deg}^{-2}$, estimated from the GOODS-ALMA survey) has
been added to each line, since $0.85\,$mJy is roughly the brightest source in
the HUDF. For each study, the contribution from stacking on undetected near-
IR-selected sources is added as a line of faded colour, arbitrarily placed at
10 $\mu$Jy. The coloured bands show two different 68 per cent confidence
estimates of the absolute value of the background.

###  5.2 The absolute value of the cosmic infrared background

Estimating the absolute level of the CIB at 1 mm is subject to larger
uncertainties than our estimate of the amount contributed by sources. The
spectrum of the CIB was measured by COBE-FIRAS (Fixsen et al., 1998), but with
fairly large systematics-dominated uncertainties. More recently Odegard et al.
(2019) used Planck in combination with FIRAS to estimate the total intensity
of the CIB in the Planck-High-Frequency Instrument channels, including at 217
and 353 GHz (the closest Planck frequencies to the mean frequency of our ALMA
image, 243 GHz). The best-fit CIB spectral shape from Fixsen et al. (1998) is
a modified blackbody function of the form

$I(\nu)=A\left(\frac{\nu}{\nu_{0}}\right)^{\beta}B_{\nu}(T_{\rm d}),$ (6)

where $A$ is a constant, $\nu_{0}$ is a fiducial frequency,
$\beta\,{=}\,0.65$, $B_{\nu}$ is the Planck function and $T_{\rm
d}\,{=}\,18.5\,$K. The uncertainties in the three fit parameters $A$, $\beta$
and $T_{\rm d}$ are significantly correlated (the correlation coefficients
reported in the paper are larger than 95 per cent), meaning that there is
effectively only one free parameter in the fit. Lacking the data needed to do
the fit ourselves, we simply fix $\beta$ and $T_{\rm d}$ (thus fixing the
shape of the CIB spectrum), but take into account the uncertainties in the
amplitude throughout our calculations.

To interpolate the CIB intensity to the frequency of our ALMA map, we first
estimate what the transmission function of our ALMA map is, and use that to
estimate the effective frequency of the image (which may be different from the
mean frequency of 243 GHz). Each ALMA observation consists of two side bands 4
GHz wide, whose central frequencies are separated by 12 GHz, and for each
observation we have already computed the central rms of the observatory-
produced MFS image ($\sigma_{i}$; see Section 3). Assuming each observation’s
individual transmission function is flat (in $\nu$) across the two sidebands,
we can calculate the mean transmission function of all of the relevant
observations used to produce our final map, weighted by each observation’s
inverse variance. To each weight we must also include a factor of the ratio of
the given observation’s beam area (set by $\theta_{{\rm maj},i}$ and
$\theta_{{\rm min},i}$) to the final image’s beam area (set by $\theta_{\rm
maj}$ and $\theta_{\rm min}$), thus the transmission function weight for
observation $i$ is

$w_{i}=\frac{1}{\sigma_{i}^{2}}\frac{\theta_{{\rm maj},i}\theta_{{\rm
min},i}}{\theta_{\rm maj}\theta_{\rm min}}.$ (7)

To obtain a transmission function, $T(\nu)$, appropriate for the average of
the entire map, we only consider the observations from surveys that uniformly
covered the HUDF, namely programmes 2012.1.00173.S (from Dunlop et al. 2017),
2015.1.00098.S (ASAGAO), 2015.1.00543.S and 2017.1.00755.S (GOODS-ALMA) and
2016.1.00324.L (ASPECS). We find that the resulting transmission function is
dominated by ASPECS and is largely flat across Band 6, with slightly more
weight towards lower frequencies.

With $T(\nu)$ in hand, we determine the effective frequency of the CIB in our
image using the modified blackbody parameters fit to the CIB spectrum from
FIRAS (Fixsen et al., 1998). Specifically, we calculate

$\nu_{\rm eff}=\frac{\int T(\nu)I(\nu)\nu d\nu}{\int T(\nu)I(\nu)d\nu},$ (8)

with $I(\nu)$ from Eq. 6 and the constants $A$ and $\nu_{0}$ cancel out in the
ratio of integrals. We find a value of $\nu_{\rm eff}\,{=}\,241.4\,$GHz, which
is close to the ALMA data combined central frequency.

We finally estimate the amplitude of the CIB at $\nu_{\rm eff}$ from Eq. 6,
using two different approaches to obtain the amplitude $A$. The first approach
is to use $A\,{=}\,(1.3\pm 0.4)\,{\times}\,10^{-5}$, which is a direct result
from Fixsen et al. (1998), assuming that the shape is fixed. The second
approach is to fit a power law between the CIB intensities from Odegard et al.
(2019) at 217 and 353 GHz (there is an exact solution because there are two
measurements and two free parameters), then use the fit to interpolate the
intensity at 241.4 GHz. To propagate uncertainties, we draw 10,000 Gaussian
random numbers from the measured amplitudes and calculate the corresponding
CIB value at 241.4 GHz.

We find absolute CIB values of $19_{-5}^{+6}\,{\rm Jy}\,{\rm deg}^{-2}$ using
the fit from Fixsen et al. (1998) and $22_{-6}^{+7}\,{\rm Jy}\,{\rm deg}^{-2}$
using the interpolation from Odegard et al. (2019), where the error bars are
the 68 per cent confidence intervals of the posterior distributions and the
central values are the means of the distributions.

We show these estimates graphically in Fig. 12. Specifically we plot with a
blue band the 68 per cent range of the absolute value of the CIB estimated
using FIRAS data alone (Fixsen et al., 1998) and with a pink band we plot the
estimate using FIRAS data in combination with Planck for foreground removal
(Odegard et al., 2019). Although the Odegard et al. (2019) results
substantially shrink some uncertainties compared to the Fixsen et al. (1998)
results, that is really only the case at shorter wavelengths; at 1.2 mm the
uncertainties are similar, but the background is actually a little higher in
the more recent paper. The difference between the two estimates indicates that
the background at these wavelengths is still quite uncertain.

In Fig. 12 we also present the contribution of sources to the CIB, by plotting
the cumulative intensity of resolved galaxies from this work as a function of
their flux density, including the estimated contribution from
$S_{1230}\,{>}\,0.85$ mJy galaxies and also the stacking-estimate contribution
from galaxies in the JADES catalogue. For reference we show the same analysis
using the results from Dunlop et al. (2017) and ASPECS (González-López et al.,
2020), where in both works an estimate of the contribution to the CIB from
stacking was performed. Our new results are similar to those from previous
studies, with our new analysis detecting the faintest sources and finding the
highest total background value in the HUDF region.

###  5.3 Consistency between resolved sources and the absolute value of the
cosmic infrared background

Clearly the absolute value measurements of the CIB are larger than what we
find by summing the flux densities of known galaxies (detected by both ALMA
and JWST). It is therefore important to consider whether or not the two kinds
of measurement are consistent.

There are fluctuations in the CIB (e.g., Planck Collaboration XXX, 2014) whose
amplitude depends on the area observed. Studies of the CIB at submm
wavelengths found that $\delta I/I\,{\simeq}\,15\,$per cent on scales around
10 arcmin (Viero et al., 2009, 2013a), which is larger than the area of our
deep ALMA map. To find a better estimate of the expected amplitude of these
fluctuations for maps the size of our HUDF image, we use the Simulated
Infrared Dusty Extragalactic Sky (SIDES) mock catalogue (Béthermin et al.,
2017; Gkogkou et al., 2023). Briefly, SIDES uses dark matter halos in a
simulated light cone to obtain clustered positions on the projected sky, then
attaches stellar masses and far-infrared SEDs to the halos using a two star-
formation-mode model of galaxy evolution. The total simulated area of the
SIDES simulation is 117 deg2, but here we only use a 1 deg2 tile from the full
simulation.

Since each simulated galaxy has an SED, we first estimate their 243-GHz flux
densities by integrating their SEDs through the transmission function derived
above. We next take 100 random patches from the 1 deg2 tile, each with the
area of our image of the HUDF (here 4.2 arcmin2). The CIB intensity of each
random patch is then just the sum of the flux densities of the galaxies in the
patch divided by the area. However, since we find no galaxies with
$S_{1230}\,{>}\,0.85$ mJy in our real ALMA map, then we subtract these from
our random patches as well; this assumes that there are no additional
fluctuations from the population of $S_{1230}\,{>}\,0.85$ mJy galaxies (i.e.
we simply add a constant for the bright part of the background). Additionally,
we neglect the contribution from galaxies at stellar masses where we were
unable to detect any statistical signal in the real data.

Following this procedure, we find that the mean CIB value within 4.2 arcmin2
patches of the SIDES simulation (after subtracting the contribution from
$S_{1230}\,{>}\,0.85$ mJy galaxies) is 11.1 Jy deg-2, comparable to the value
of $(8.6\,{\pm}\,0.5)\,{\rm Jy}\,{\rm deg}^{-2}$ that we have measured. The
standard deviation of the SIDES simulation patches is 1.9 Jy deg-2, thus
fluctuations are $\delta I/I\,{=}\,17\,$per cent. We must therefore take this
into account when comparing the mean CIB value from FIRAS averaged over nearly
the whole sky to the single 4.2 arcmin2 patch we have observed.

We now estimate the probability of measuring a CIB value of
$(10.0\,{\pm}\,0.5)\,{\rm Jy}\,{\rm deg}^{-2}$ assuming that the true value is
either $19_{-5}^{+6}\,{\rm Jy}\,{\rm deg}^{-2}$ (Fixsen et al., 1998) or
$23_{-8}^{+6}\,{\rm Jy}\,{\rm deg}^{-2}$ (Odegard et al., 2019) and that
fluctuations are 17 per cent. We take our 10,000 CIB absolute value
realizations discussed above and additionally draw 10,000 Gaussian-distributed
values for a fluctuation amplitude, then multiply these together. We draw
10,000 Gaussian-distributed numbers for our measured CIB value, and take the
difference between these and the possible absolute values. Finally, our
statistic is the fraction of the area of the resulting posterior distribution
that is less than zero, which can be interpreted as the probability of
obtaining our actual measured CIB value or less while taking into account both
the measurement uncertainties and intrinsic CIB fluctuations. We find that the
probability of measuring a CIB value of $(10.0\,{\pm}\,0.5)\,{\rm Jy}\,{\rm
deg}^{-2}$, given the absolute CIB measurement from Fixsen et al. (1998) is
8.9 per cent, or 5.3 per cent given the absolute CIB measurement from Odegard
et al. (2019).

Assuming that the absolute value measurements from FIRAS are correct and that
our measurement is also correct, we can calculate the required level of
statistical excursion (in units of $\sigma\,{=}\,0.17$) corresponding to the
HUDF. Using the same 10,000 random values for the absolute FIRAS values and
the measured values, we find that the CIB fluctuation at the position of the
HUDF must be $-2.4_{-1.4}^{+0.4}\sigma$ using the CIB value from Fixsen et al.
(1998), or $-2.9_{-1.2}^{+0.3}\sigma$ from Odegard et al. (2019) (here the
central values are the means of the posterior distributions and the error bars
are 68 per cent confidence intervals). What this means is that the variance in
HUDF fields is large enough that our results can explain the whole of the CIB
provided that the HUDF happens to be a relatively mild (${\simeq}\,2\,\sigma$)
underdense direction on the sky.666The HUDF was not selected entirely
randomly; however, there is no particular reason to believe that the criteria
used for its selection would make it likely to have a lower than average
background (Beckwith et al., 2006).

As a final check, we use the simulated catalogue of galaxies from SIDES to
estimate the fraction of the CIB emitted by galaxies with stellar masses
between 107.4 M⊙ and 1011.4 M⊙ and with redshifts between 0 and 8 (namely the
parameter space over which we stacked on undetected JADES galaxies). For each
of the 100 random patches described above, we also calculate the sum of the
flux densities from all galaxies within our stacking range divided by the sum
of all of the galaxies in the HUDF-sized region. We find that the average
ratio is 97 per cent, with a standard deviation of 1 per cent. Thus if we
accept that SIDES provides a reasonable model for counts at these wavelengths,
we do not expect that we are missing a significant contribution to our ALMA
measurement of the CIB from even fainter and lower stellar mass galaxies.

Turning back to the data, if we include sources detected by ALMA, the
contribution to the CIB from 1010.4–1011.4 M⊙ galaxies is
$(3.5\,{\pm}\,0.1)\,{\rm Jy}\,{\rm deg}^{-2}$, from 109.4–1010.4 M⊙ galaxies
is $(3.2\,{\pm}\,0.1)\,{\rm Jy}\,{\rm deg}^{-2}$, from 108.4–109.4 M⊙ galaxies
is $(1.1\,{\pm}\,0.2)\,{\rm Jy}\,{\rm deg}^{-2}$ and from 107.4–108.4 M⊙
galaxies is $(0.5\,{\pm}\,0.4)\,{\rm Jy}\,{\rm deg}^{-2}$. If we stack on the
next lowest mass bin (106.4–107.4 M⊙) over the full redshift range we obtain a
CIB contribution of $(-0.5\,{\pm}\,0.5)\,{\rm Jy}\,{\rm deg}^{-2}$. This is
consistent with the CIB having essentially converged over the stellar mass
range that we have probed.

A potential issue here is the completeness of the JADES catalogue within our
stacking range. For our 1000 random SIDES patches we also keep track of the
total number of galaxies with stellar masses and redshifts within our stacking
region. We find that the average number of galaxies is 1890 (with a standard
deviation of 150). In the JADES catalogue there are 1856 galaxies in our ALMA
image with stellar masses and redshifts within our stacking range, which is
consistent with the total number of galaxies expected from the SIDES
simulation. For reference, there are 1561 galaxies from the CANDELS catalogue
that are in our ALMA image footprint within the same stacking range. Adding
JADES has helped us to find more of the CIB within our ALMA map, and it seems
that going even deeper in the optical/near-IR will not add significantly to
the source-derived CIB estimate.

What these numbers indicate is that our estimate of the 1.23-mm CIB (from
individually-detected galaxies, together with statistical stacking results)
appears to contain essentially all of the possible galaxies that would
contribute to the CIB, and that it is genuinely lower than what the mean value
is estimated to be. However, the chances that our small patch of the
extragalactic sky falls on a negative fluctuation of the CIB are not small
enough to rule out the hypothesis that we have indeed recovered essentially
the entire CIB from these known galaxies.

##  6 Improvement over previous studies

The deepest previously-published ALMA survey of the HUDF at 1 mm is ASPECS and
so a question of interest is how much of an improvement have we achieved by
including all of the additional data. Qualitatively, looking at Table 1 we can
see that ASPECS is by far the deepest of the individual surveys. The maps from
Dunlop et al. (2017), ASAGAO and GOODS-ALMA overlap with the entire ASPECS map
and so by combining them with ASPECS the result must be deeper. Additionally,
by including more individual pointings, we have been able to make some regions
even deeper still.

To quantify the difference between ASPECS and our combined image, we
downloaded the ASPECS map produced using CASA in the MFS mode and made public
by the ASPECS team,777https://almascience.nrao.edu/alma-data/lp/ASPECS then
ran the primary-beam-corrected image through our algorithm for generating the
noise map (see Section 2.3). The pixel size of the ASPECS map is 0.2 arcsec,
the same as our $uv$-combined map, and the beamsizes are very similar
($1.49\,{\rm arcsec}\,{\times}\,1.07\,{\rm arcsec}$ versus $1.53\,{\rm
arcsec}\,{\times}\,1.08\,{\rm arcsec}$), so we simply computed the ratio of
the two noise maps to assess the amount by which the noise improves with the
additional data. Unsurprisingly we find that across the entire ASPECS region
the noise (meaning here the rms after masking sources) in our map is smaller,
ranging from about 5 per cent smaller in the central region to about 50 per
cent smaller near the edges and around the deepest individual pointings. The
ASPECS map was made going out to a primary beam value of 0.1 and covers a
total area of 4.2 arcmin2, 3.3 arcmin2 of which has a noise level less than 35
$\mu$Jy beam-1, whereas our map was made going out to a primary beam value of
0.2 and all 4.2 arcmin2 has a noise level less than 35 $\mu$Jy beam-1. Thus by
combining the different data sets we not only reduce the noise across the
majority of the map by 5 per cent, but we also expand the area where the map
is at its most sensitive.

González-López et al. (2020) searched the ASPECS map for sources down to a
fidelity of 0.5, and their lowest significance source had a S/N of 3.3. Their
catalogue contains a total of 35 sources (four of which are not detected in
our map), whereas we find 45 sources; it is worth noting that all 45 of our
sources are contained within the area mapped by ASPECS. Thus the extra depth
we have been able to achieve by combining archival ALMA data with ASPECS has
led to a significant increase in detected sources.

In addition to detecting more sources, we are able to recover more of the CIB
through stacking thanks to our deeper JADES catalogue. If we stack on the 1661
galaxies from the CANDELS catalogue with stellar masses between 107.4 M⊙ and
1011.4 M⊙ and with redshifts between 0 and 8 (after masking detected ALMA
sources) we obtain $(1.5\,{\pm}\,0.4)\,{\rm Jy}\,{\rm deg}^{-2}$, compared to
$(2.4\,{\pm}\,0.5)\,{\rm Jy}\,{\rm deg}^{-2}$ using the JADES catalogue.

Now turning to the wider HUDF region, Hatsukade et al. (2018) also combined
the ASAGAO survey with the survey from Dunlop et al. (2017) and the first
GOODS-ALMA survey (Franco et al., 2018); their final map has a mean rms of
about $75\,\mu{\rm Jy}\,{\rm beam}^{-1}$, with a beamsize of $0.59\,{\rm
arcsec}\,{\times}\,0.53\,{\rm arcsec}$ and a pixel size of 0.1 arcsec after
applying a 250 k$\lambda$ taper. Our shallow and wide combined map (presented
in Appendix A) contains additional GOODS-ALMA data that were not available
when the ASAGAO map was constructed, as well as more individual pointings,
while maintaining a similar beamsize ($0.87\,{\rm
arcsec}\,{\times}\,0.64\,{\rm arcsec}$) and pixel size (0.12 arcsec). We thus
also downloaded the ASAGAO 250 k$\lambda$-tapered
map888https://sites.google.com/view/asagao26/alma-data?authuser=0 to
quantitatively check the improvement. After running the ASAGAO map through our
noise algorithm, we find that the new GOODS-ALMA data reduces the noise by
10–20 per cent throughout, while some individual pointings go about 50 per
cent deeper.

##  7 Conclusions

We have produced a series of 1-mm maps of the HUDF by combining all of the
previously-published survey data in the $uv$ plane in various ways, reducing
the noise compared to previous studies. We specifically constructed a deep map
covering 4.2 arcmin2 and a shallower map covering 25.4 arcmin2. Our deep map
has a pixel rms that ranges from 5 to 50 per cent lower than in the best
previous study, with an area of about 1.5 arcmin reaching below $9\,\mu{\rm
Jy}\,{\rm beam}^{-1}$ and a minimum of about $4.6\,\mu{\rm Jy}\,{\rm
beam}^{-1}$ reached in some regions. Our shallow map has a pixel rms that
ranges from 10 to 50 per cent lower than in the best previous wider and
shallower study. We make all of our maps publicly
available999https://doi.org/10.5683/SP3/YWBVWH

We searched our deep map for sources down to a signal-to-noise threshold of
3.6, finding a total of 45 peaks in the S/N map, 13 of which are new. Nearly
all (38/45) of these ALMA sources have near-IR counterparts detected by JWST
and HST. We additionally find 27 sources in our wider map, nine of which are
new. The JADES data enable stellar masses and photometric redshifts to be
estimated for the ALMA source counterparts, and we find that they are all
relatively high $M_{\ast}$ galaxies. Compared with ALMA-undetected galaxies at
similar $M_{\ast}$, the ALMA-detected galaxies typically have redder JWST
colours. With our larger sample of mm-selected sources in the HUDF, other
studies investigating the statistical properties of the faintest star-forming
galaxies could be carried out. For example, Aravena et al. (2020) studied the
SFR-stellar mass relation of ASPECS galaxies and Boogaard et al. (2023) looked
at the morphologies of ASPECS galaxies in JWST’s MIRI filters; expanding to a
larger sample size and using improved SED fits from JWST photometry could lead
to more robust conclusions.

Since the vast majority of near-IR-selected galaxies are not directly detected
in our 1-mm map, we performed a stacking analysis on their positions. We found
significant average signals for all galaxies in the range $z\,{=}\,0$ to 3 and
with stellar masses between $10^{9.4}\,{\rm M}_{\odot}$ and $10^{10.4}\,{\rm
M}_{\odot}$, as well as a roughly $3\,\sigma$ signal from $10^{8.4}\,{\rm
M}_{\odot}$ to $10^{9.4}\,{\rm M}_{\odot}$ stellar mass galaxies between
$z\,{=}\,0$ and 2. The evidence is that $M_{\ast}\,{\sim}\,10^{10}\,{\rm
M}_{\odot}$ ALMA-undetected galaxies have a 1-mm flux density around
$15\,\mu$Jy and would be individually detected in even deeper integrations.
The stacking results also show the value of our new data products for
performing similar statistical analyses on sets of galaxies detected in other
wavebands.

We used our galaxy detections, as well as our stacking analysis, to estimate
the level of the CIB at 1.23 mm that we have resolved. We thus account for a
background level of $(10.0\,{\pm}\,0.5)\,{\rm Jy}\,{\rm deg}^{-2}$, with an
expectation that even fainter galaxies will hardly change this number. There
is still a large uncertainty in the total background estimate, and we also
stress that the variance expected in a region as small as the HUDF is around
20 per cent. We do not recover all of the background, and there are a number
of possible resolutions: (1) HUDF is a 2–$3\,\sigma$ negative fluctuation in
the CIB; (2) the absolute level of the CIB has been overestimated; (3) there
is a new population of galaxies that contribute at 1 mm, but are not detected
by JWST; or (4) a genuinely diffuse component of the background exists. Given
that the simple explanation (1) cannot be excluded, then the suggestion is
that deep ALMA+JWST observations may account for all of the mm-wave
background.

The HUDF is probably the best-studied extragalactic field, and it is crucial
to probe this field at longer (${>}\,500\,\mu$m) wavelengths in order to
understand how the earliest galaxies began forming their stars, complementing
the data obtained by HST, JWST and other telescopes at shorter wavelengths.
The 4.2 arcmin2 map presented here is the deepest such image of the HUDF
available to-date, yet it is clear that there are still many more galaxies
left to detect. Deeper ALMA observations of the HUDF will inevitably uncover
these galaxies, providing a more complete understanding of galaxy evolution at
early times.

## Acknowledgements

This research used the Canadian Advanced Network For Astronomy Research
(CANFAR) operated in partnership by the Canadian Astronomy Data Centre and The
Digital Research Alliance of Canada, with support from the National Research
Council of Canada, the Canadian Space Agency, CANARIE and the Canadian
Foundation for Innovation. RH and DS acknowledge the support of the Natural
Sciences and Engineering Research Council of Canada (NSERC). This paper makes
use of the ALMA data ADS/JAO.ALMA#2012.1.00173.S, 2013.1.00718.S,
2013.1.01271.S, 2015.1.00098.S, 2015.1.00543.S, 2015.1.00664.S,
2015.1.00821.S, 2015.1.00870.S, 2015.1.01096.S, 2015.1.01379.S,
2015.1.01447.S, 2015.A.00009.S, 2016.1.00324.L, 2016.1.00721.S,
2016.1.00967.S, 2017.1.00190.S, 2017.1.00755.S, 2018.1.00567.S and
2018.1.01044.S. ALMA is a partnership of ESO (representing its member states),
NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA
(Taiwan) and KASI (Republic of Korea), in cooperation with the Republic of
Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. This
research made use of Photutils, an Astropy package for detection and
photometry of astronomical sources. This work is based in part on observations
taken by the CANDELS Multi-Cycle Treasury Program with the NASA/ESA HST, which
is operated by the Association of Universities for Research in Astronomy,
Inc., under NASA contract NAS5-26555. This work is based in part on
observations made with the NASA/ESA/CSA James Webb Space Telescope.

## Data Availability

All of the data products described in this paper are publicly available at
https://doi.org/10.5683/SP3/YWBVWH. The raw ALMA data used to produce these
data products are publicly available at the ALMA archive. The JWST images from
the JADES Collaboration used in this paper are also publicly available, and
the catalogue derived using JADES images are described and made available in a
separate paper.

## References

  * Aravena et al. (2016) Aravena M., et al., 2016, ApJ, 833, 68
  * Aravena et al. (2020) Aravena M., et al., 2020, ApJ, 901, 79
  * Arnouts & Ilbert (2011) Arnouts S., Ilbert O., 2011, LePHARE: Photometric Analysis for Redshift Estimate, Astrophysics Source Code Library, record ascl:1108.009 (ascl:1108.009)
  * Beckwith et al. (2006) Beckwith S. V. W., et al., 2006, AJ, 132, 1729
  * Bertin & Arnouts (1996) Bertin E., Arnouts S., 1996, A&AS, 117, 393
  * Béthermin et al. (2017) Béthermin M., et al., 2017, A&A, 607, A89
  * Boogaard et al. (2023) Boogaard L. A., et al., 2023, arXiv e-prints, p. arXiv:2308.16895
  * Borys et al. (2003) Borys C., Chapman S., Halpern M., Scott D., 2003, MNRAS, 344, 385
  * Bouwens et al. (2011) Bouwens R. J., et al., 2011, ApJ, 737, 90
  * Bradley et al. (2022) Bradley L., et al., 2022, astropy/photutils: 1.5.0, https://doi.org/10.5281/zenodo.6825092
  * Bruzual & Charlot (2003) Bruzual G., Charlot S., 2003, MNRAS, 344, 1000
  * Calzetti et al. (2000) Calzetti D., Armus L., Bohlin R. C., Kinney A. L., Koornneef J., Storchi-Bergmann T., 2000, ApJ, 533, 682
  * Casey et al. (2015) Casey C. M., et al., 2015, ApJ, 808, L33
  * Chabrier (2003) Chabrier G., 2003, PASP, 115, 763
  * Chen et al. (2023) Chen J., et al., 2023, MNRAS, 518, 1378
  * Condon (1997) Condon J. J., 1997, PASP, 109, 166
  * Condon et al. (1998) Condon J. J., Cotton W. D., Greisen E. W., Yin Q. F., Perley R. A., Taylor G. B., Broderick J. J., 1998, AJ, 115, 1693
  * Cowie et al. (2018) Cowie L. L., González-López J., Barger A. J., Bauer F. E., Hsu L. Y., Wang W. H., 2018, ApJ, 865, 106
  * Dunlop et al. (2017) Dunlop J. S., et al., 2017, MNRAS, 466, 861
  * Eisenstein et al. (2023) Eisenstein D. J., et al., 2023, arXiv e-prints, p. arXiv:2306.02465
  * Fixsen et al. (1998) Fixsen D. J., Dwek E., Mather J. C., Bennett C. L., Shafer R. A., 1998, ApJ, 508, 123
  * Franco et al. (2018) Franco M., et al., 2018, A&A, 620, A152
  * Fujimoto et al. (2017) Fujimoto S., Ouchi M., Shibuya T., Nagai H., 2017, ApJ, 850, 83
  * Gkogkou et al. (2023) Gkogkou A., et al., 2023, A&A, 670, A16
  * Gómez-Guijarro et al. (2022) Gómez-Guijarro C., et al., 2022, A&A, 658, A43
  * González-López et al. (2020) González-López J., et al., 2020, ApJ, 897, 91
  * Guo et al. (2013) Guo Y., et al., 2013, ApJS, 207, 24
  * Hatsukade et al. (2018) Hatsukade B., et al., 2018, PASJ, 70, 105
  * Hill et al. (2018) Hill R., Masui K. W., Scott D., 2018, Applied Spectroscopy, 72, 663
  * Illingworth et al. (2013) Illingworth G. D., et al., 2013, ApJS, 209, 6
  * Ivison et al. (2007) Ivison R. J., et al., 2007, MNRAS, 380, 199
  * Koekemoer et al. (2013) Koekemoer A. M., et al., 2013, ApJS, 209, 3
  * Koprowski et al. (2017) Koprowski M. P., Dunlop J. S., Michałowski M. J., Coppin K. E. K., Geach J. E., McLure R. J., Scott D., van der Werf P. P., 2017, MNRAS, 471, 4155
  * Major et al. (2019) Major B., Kavelaars J., Fabbro S., Durand D., Jeeves H., 2019, in Teuben P. J., Pound M. W., Thomas B. A., Warner E. M., eds, Astronomical Society of the Pacific Conference Series Vol. 523, Astronomical Data Analysis Software and Systems XXVII. p. 277
  * Marsden et al. (2009) Marsden G., et al., 2009, ApJ, 707, 1729
  * Massardi et al. (2021) Massardi M., et al., 2021, PASP, 133, 085001
  * McLure et al. (2018) McLure R. J., et al., 2018, MNRAS, 476, 3991
  * McMullin et al. (2007) McMullin J. P., Waters B., Schiebel D., Young W., Golap K., 2007, in Shaw R. A., Hill F., Bell D. J., eds, Astronomical Society of the Pacific Conference Series Vol. 376, Astronomical Data Analysis Software and Systems XVI. p. 127
  * Merlin et al. (2015) Merlin E., et al., 2015, A&A, 582, A15
  * Nonino et al. (2009) Nonino M., et al., 2009, ApJS, 183, 244
  * Odegard et al. (2019) Odegard N., Weiland J. L., Fixsen D. J., Chuss D. T., Dwek E., Kogut A., Switzer E. R., 2019, ApJ, 877, 40
  * Penner et al. (2011) Penner K., et al., 2011, MNRAS, 410, 2749
  * Planck Collaboration XXX (2014) Planck Collaboration XXX 2014, A&A, 571, A30
  * Pope et al. (2005) Pope A., Borys C., Scott D., Conselice C., Dickinson M., Mobasher B., 2005, MNRAS, 358, 149
  * Rieke & the JADES Collaboration (2023) Rieke M., the JADES Collaboration 2023, arXiv e-prints, p. arXiv:2306.02466
  * Robitaille et al. (2020) Robitaille T., Deil C., Ginsburg A., 2020, reproject: Python-based astronomical image reprojection, Astrophysics Source Code Library, record ascl:2011.023 (ascl:2011.023)
  * Viero et al. (2009) Viero M. P., et al., 2009, ApJ, 707, 1766
  * Viero et al. (2013a) Viero M. P., et al., 2013a, ApJ, 772, 77
  * Viero et al. (2013b) Viero M. P., et al., 2013b, ApJ, 779, 32
  * Whitaker et al. (2019) Whitaker K. E., et al., 2019, ApJS, 244, 16
  * Williams et al. (2023) Williams C. C., et al., 2023, arXiv e-prints, p. arXiv:2301.09780

##  Appendix A Sources in the wider region

Table 4: Archival ALMA projects downloaded from the ALMA archive and used to
create the larger and shallower combined 1-mm image.

Project ID | Target name(s) | Frequency range | Map rmsa | Synthesized beamb | Sky coveragec
---|---|---|---|---|---
|  | [GHz] | [$\mu$Jy beam-1] | [arcsec$\,{\times}\,$arcsec] | [arcmin2]
2015.1.00098.Sd | HUDF-JVLA-ALMA | 244.3–271.8 | 57 | 0.20$\,{\times}\,$0.16 | 34.7
2015.1.00543.Sd | GOODS-S | 255.1–274.7 | 130 | 0.26$\,{\times}\,$0.22 | 54.2
2015.1.00664.S | KMOS3DGS4-27882 | 255.2–274.8 | 63 | 0.15$\,{\times}\,$0.14 | 0.23
2015.1.00821.S | z7_GSD_3811 | 218.2–236.6 | 23 | 0.70$\,{\times}\,$0.51 | 0.32
2015.1.00870.S | 19842 | 223.2–242.8 | 42 | 0.58$\,{\times}\,$0.49 | 0.30
2015.1.01379.S | GMASS_0953 | 212.3–216.8 | 26 | 0.66$\,{\times}\,$0.54 | 0.36
2016.1.00721.S | z6_GSD_1388 | 250.9–269.3 | 23 | 0.78$\,{\times}\,$0.53 | 0.24
2016.1.00967.S | ZFOURGE_CDSF_4409 | 243.8–262.1 | 28 | 1.24$\,{\times}\,$0.84 | 0.26
2017.1.00755.Sd | GOODS-S | 255.1–274.9 | 120 | 1.31$\,{\times}\,$0.82 | 54.1
2017.1.00190.S | z7_GSD_3811 | 218.4–236.9 | 12 | 0.61$\,{\times}\,$0.49 | 0.32
2018.1.00567.S | All ASAGAO pointings | 244.2–262.9 | 29 | 0.65$\,{\times}\,$0.47 | 5.15
2018.1.01044.S | cdfs_535, cdfs_769 | 223.1–242.9 | 27 | 0.63$\,{\times}\,$0.51 | 0.61

  * •

a Approximate map rms estimated using the observatory MFS data products. For
programmes with a single tuning, this is the rms of the primary-beam-
uncorrected map after masking all known sources. For programmes with multiple
tunings, we follow the same procedure and then estimate the rms of the
weighted mean of the images as $1/\big{(}\Sigma_{i}1/\sigma_{i}^{2}\big{)}$.
The actual rms from combining images at different tunings in $uv$ space using
the MFS mode will generally be less than this estimate.

  * •

b Maximum synthesized beam across all frequencies of a given data set. The
actual synthesized beam from combining data taken in different tunings using
the MFS mode may differ slightly from this estimate.

  * •

c The sky coverage overlapping with our data selection criteria, which may be
less than the total sky coverage of the given programme.

  * •

d Programmes used to make a mosaic with uniform noise.

Figure 13: Left: Signal map of the shallower, larger region of the HUDF after
combining all of the available archival ALMA Band-6 data given in Table 4.
This region is defined by the contour where the primary beam reaches 0.4,
covering a total of 25.4 arcmin2. Middle: Corresponding noise map. Right: S/N
map, made by dividing the signal map by the noise map. In all panels, the cyan
contour shows the footprint of our combined deeper central map, which we
exclude from our source extraction procedure in the shallower region. Figure
14: S/N map of the larger and shallower combined ALMA data, covering 25.4
arcmin2, with peaks ${>}\,4.5\,\sigma$ (where the fidelity is 0.9) indicated
as green boxes. Galaxies from the ASAGAO survey (Hatsukade et al., 2018) are
shown as red circles and galaxies found in the GOODS-ALMA survey (Gómez-
Guijarro et al., 2022) are shown as brown circles. Galaxies with counterparts
in both the CANDELS catalogue and the JADES catalogue are indicated with a
blue cross, galaxies with only a CANDELS counterpart are indicated with a gold
cross and galaxies with only a JADES counterpart are indicated with a red
cross. The dashed cyan contour shows the footprint of our combined deeper
central map, which we exclude from our source extraction procedure in this
appendix, while the gold contour shows the footprint of the deepest HST data
from the CANDELS survey and the red contour shows the footprint of the JADES
survey. Figure 15: Comparison of the flux densities measured in the outer
region of the map versus flux densities from ASAGAO (Hatsukade et al., 2018)
and GOODS-ALMA (Gómez-Guijarro et al., 2022).

Here we present a shallower but larger ALMA 1-mm map, in an area defined by
the ASAGAO footprint, minus the footprint of the deeper inner region. We will
also give a catalogue of sources found in this wider map. The total area of
the shallower map is 25.4 arcmin2, so after subtracting the area of the inner
region we end up with an effective area of 21.2 arcmin2. The typical pixel rms
is $60\,\mu{\rm Jy}\,{\rm beam}^{-1}$ outside of the individual pointings and
$20\,\mu{\rm Jy}\,{\rm beam}^{-1}$ within the individual pointings.

Table 4 outlines the programmes used to produce the combined image, Fig. 13
shows the resulting signal, noise and S/N maps, Table 5 provides positions,
flux densities, redshifts and stellar masses for all of the sources found
above a threshold of $4.5\,\sigma$ (where the fidelity reaches a value of
0.9), ignoring the deeper central region where we can find fainter sources;
for reference, the most significant negative peak found in this search has a
(nagative) S/N of 5.1. Finally, Table 6 provides IDs matched to other
catalogues. One ALMA source (ID 51) is not included in the JADES catalogue
(likely due to blending with another nearby galaxy), so we list the ID from
McLeod et al. (in prep.) instead.

Flux densities were measured using the same method described in Section 4.1,
with one adjustment. We found that a 2D Gaussian was a poor fit to high S/N
sources (presumably mostly because of the combination of data sets with
different resolution), which resulted in flux densities a factor of about 2
larger than reported in previous studies. We therefore set a threshold of 1000
$\mu$Jy above which we do not fit a 2D Gaussian, and instead just take the
peak pixel value. We found that this resulted in flux densities consistent
with previous studies. Lastly, our source ALMA-HUDF-69 has a counterpart ID 44
in the ASAGAO catalogue, where it is reported that the peak pixel flux density
is nearly 6 times fainter than the integrated flux density. The JWST/HST
counterpart to this galaxy is a large resolved spiral at redshift 1.69, so it
is likely that this galaxy is also highly resolved in our ALMA image. We find
that a 2D Gaussian is a poor fit to the ALMA data, so we simply use an
aperture with a radius of 2 arcsec to measure its flux density. As an
additional test, we extracted peaks from the central region of the shallow map
and measured their flux densities in order to compare them with our deep map;
we found good agreement between the two measurements for all ${>}\,4.5\sigma$
sources in the shallow map.

In Fig. 14 we show the positions of our sources, along with sources found in
the ASAGAO (Hatsukade et al., 2018) and GOODS-ALMA (Gómez-Guijarro et al.,
2022) surveys. Blue crosses indicate sources with a CANDELS and JADES match,
gold crosses indicate sources with a CANDELS-only match (these galaxies all
lie outside of the JADES footprint) and red crosses indicate galaxies with a
JADES-only match. Figure 15 compares our photometry to the photometry of the
ASAGAO and GOODS-ALMA surveys for all overlapping sources. In Appendix B we
show cutouts of these ALMA sources overlaid on JWST F356W imaging.

We find 27 galaxies in our defined region, nine of which are not in a
previously-published catalogue. Of these nine, five have a JADES/CANDELS
counterpart, while the remaining four do not; however, three of these four
galaxies lie outside of both the deep HST footprint and the JADES footprint.
It may also be worth noting that the fourth brightest source (ALMA-HUDF-51) in
this supplementary catalogue lies just beyond the boundary of our deep map.

We repeated the analyses of the sources in the wide map that we performed for
the deep map. In particular we looked at: the redshift, magnitude and colour
distributions of the previously known versus new sources, as in Fig. 7; the
stellar mass versus redshift distribution, as in Fig. 8; and the magnitude and
colour distributions of the ALMA galaxies versus ALMA-undetected galaxies of
similar stellar mass, as in Fig. 9. In each case the results were similar to
those from the deep map, but with fewer sources.

It would also be possible to carry out a stacking analysis in the wider
region. However, for a fixed number density of optical/near-IR galaxies, the
stacking uncertainty in a bin, $\delta\hat{S}_{\alpha}$, depends on the survey
area and depth of the 1-mm map in a way that essentially scales like the
square root of the total observing time. We can see that in the following way.
Let us assume that we have a catalogue with a fixed number density of sources,
$n_{\rm s}$, and we are stacking on two maps, one with area $A_{1}$ and noise
per pixel (or beamsized-pixel, say) $\sigma_{1}$, and the other with area
$A_{2}$ and noise per pixel $\sigma_{2}$. The total numbers of sources to
stack on in the two maps are then $N_{1}=n_{\rm s}A_{1}$ and $N_{2}=n_{\rm
s}A_{2}$. The error on the stack will be the noise in a pixel reduced by the
square root of the number of samples, hence the ratio if errors between the
two maps will be

$\frac{\delta\hat{S}_{1}}{\delta\hat{S}_{2}}=\frac{\sigma_{1}}{\sigma_{2}}\left(\frac{N_{2}}{N_{1}}\right)^{1/2}=\frac{\sigma_{1}}{\sigma_{2}}\left(\frac{n_{\rm
s}A_{2}}{n_{\rm
s}A_{1}}\right)^{1/2}=\frac{\sigma_{1}}{\sigma_{2}}\left(\frac{A_{2}}{A_{1}}\right)^{1/2}.$
(9)

But the time taken to map an area to fixed depth is proportional to the area,
while the time to map a fixed area to a given depth is proportional to the
inverse square of the depth. In other words the integration times are
$t_{1}\propto A_{1}/\sigma_{1}^{2}$ and $t_{2}\propto A_{2}/\sigma_{1}^{2}$,
which means that the ratio in Eq. 9 is
$\delta\hat{S}_{1}/\delta\hat{S}_{2}=(t_{2}/t_{1})^{1/2}$. For the case at
hand, since so much more integration time went into the deep map, then
stacking results from the wider map add will have much larger uncertainties
and provide very little additional information.

Table 5: Positions and flux densities ($S_{1230}$) for the 27 sources found in the ALMA 1-mm combined image. Here flux densities are measured by fitting Gaussian profiles to the sources in order to estimate the number of beams per source, then scaling the peak pixel value by the number of beams. For sources consistent with 1 beam or less and for sources with peak pixels ${>}\,1000\,\mu$Jy we set the flux density to be the peak pixel value. For the extended galaxy ALMA-HUDF-69 (marked with a dagger), we use an aperture with a radius of 2 arcsec to calculate a flux density. For galaxies with matches in the JADES catalogue and sufficient photometry to fit SEDs we provide stellar masses and photometric redshifts from McLeod et al. (in prep.). Galaxies with spectroscopic redshifts from Rieke & the JADES Collaboration (2023) are indicated by asterisks, which we provide instead of photometric redshifts. The median photometric redshift uncertainty of our ALMA sample is $\Delta z/(1\,{+}\,z)\,{\simeq}\,0.03$, while the relative uncertainties in stellar masses are small. Name | RA Dec [J2000] |  $S_{1230}$ [$\mu$Jy] |  $z$0 | $\log\left(M_{\ast}/\mathrm{M}_{\odot}\right)$
---|---|---|---|---
ALMA-HUDF-46 | 3:32:28.51 $-$27:46:58.3 | 1976$\pm$260 | …0 | …
ALMA-HUDF-47 | 3:32:35.73 $-$27:49:16.2 | 1870$\pm$750 | 2.410∗ | 11.1
ALMA-HUDF-48 | 3:32:34.27 $-$27:49:40.4 | 1863$\pm$820 | 3.580∗ | 10.2
ALMA-HUDF-49 | 3:32:44.04 $-$27:46:35.9 | 1146$\pm$710 | 2.850∗ | 10.8
ALMA-HUDF-50 | 3:32:28.59 $-$27:48:50.5 |  0939$\pm$157 | 3.594∗ | 10.4
ALMA-HUDF-51 | 3:32:32.91 $-$27:45:41.0 |  0890$\pm$260 | 1.160∗ |  09.6
ALMA-HUDF-52 | 3:32:34.01 $-$27:49:00.0 |  0889$\pm$140 | 2.630∗ | 10.9
ALMA-HUDF-53 | 3:32:28.78 $-$27:44:35.3 |  0879$\pm$750 | …0 | …
ALMA-HUDF-54 | 3:32:33.25 $-$27:49:16.5 |  0828$\pm$173 | 3.560∗ |  09.5
ALMA-HUDF-55 | 3:32:48.00 $-$27:46:27.0 |  0811$\pm$181 | 2.03 | 10.2
ALMA-HUDF-56 | 3:32:47.60 $-$27:44:52.5 |  0768$\pm$650 | 1.70 | 10.6
ALMA-HUDF-57 | 3:32:28.91 $-$27:44:31.6 |  0719$\pm$440 | …0 | …
ALMA-HUDF-58 | 3:32:49.46 $-$27:49:08.9 |  0717$\pm$390 | …0 | …
ALMA-HUDF-59 | 3:32:30.36 $-$27:45:01.4 |  0680$\pm$144 | …0 | …
ALMA-HUDF-60 | 3:32:31.48 $-$27:46:23.4 |  0671$\pm$300 | 2.225∗ | 10.9
ALMA-HUDF-61 | 3:32:44.60 $-$27:48:36.2 |  0642$\pm$115 | 2.593∗ | 10.7
ALMA-HUDF-62 | 3:32:29.25 $-$27:45:10.0 |  0568$\pm$570 | …0 | …
ALMA-HUDF-63 | 3:32:45.19 $-$27:48:06.9 |  0553$\pm$830 | 3.940∗ |  09.5
ALMA-HUDF-64 | 3:32:42.15 $-$27:48:47.2 |  0493$\pm$112 | …0 | …
ALMA-HUDF-65 | 3:32:47.17 $-$27:45:25.5 |  0478$\pm$520 | 3.620∗ | 10.2
ALMA-HUDF-66 | 3:32:40.76 $-$27:49:26.4 |  0441$\pm$910 | 2.170∗ | 10.6
ALMA-HUDF-67 | 3:32:28.83 $-$27:48:29.9 | 1255$\pm$382 | 1.687∗ | 11.0
ALMA-HUDF-68 | 3:32:38.74 $-$27:48:40.1 |  0389$\pm$630 | 2.640∗ | 10.9
ALMA-HUDF-69 | 3:32:43.68 $-$27:48:51.0 |  0303$\pm$600 | 2.400∗ | 10.3
ALMA-HUDF-70 | 3:32:48.60 $-$27:48:10.3 |  0300$\pm$640 | …0 | …
ALMA-HUDF-71 | 3:32:30.29 $-$27:45:24.9 |  0252$\pm$550 | …0 | …
ALMA-HUDF-72 | 3:32:41.76 $-$27:48:25.1 |  0179$\pm$440 | 1.570∗ | 10.1
Table 6: Cross-matched IDs between our catalogue of sources and previous studies. Near-infrared matches are from JADES (Rieke & the JADES Collaboration, 2023) and CANDELS (Guo et al., 2013). Matches with other ALMA Band-6 surveys are given from ASAGAO (Hatsukade et al., 2018) and GOODS-ALMA (Gómez-Guijarro et al., 2022). One ALMA source (ID 51) does not have a JADES ID (likely due to blending with another nearby galaxy), so we provide the ID from from McLeod et al. (in prep.) instead. Name | JADES ID | CANDELS ID | ASAGAO ID | GOODS-ALMA ID
---|---|---|---|---
ALMA-HUDF-46 | … | J033228.50$-$274658.1 | 2 | A2GS1
ALMA-HUDF-47 | 196290 | J033235.71$-$274916.0 | 3 | A2GS3
ALMA-HUDF-48 | 193893 | … | … | A2GS2
ALMA-HUDF-49 | 209026 | J033244.02$-$274635.7 | 1 | A2GS21
ALMA-HUDF-50 | 198459 | … | 33 | …
ALMA-HUDF-51 | 38804† | … | 7 | A2GS28
ALMA-HUDF-52 | 197581 | … | 31 | …
ALMA-HUDF-53 | … | J033228.77$-$274434.9 | … | A2GS37
ALMA-HUDF-54 | 96357 | J033233.24$-$274916.0 | … | …
ALMA-HUDF-55 | 209737 | J033247.99$-$274626.6 | … | …
ALMA-HUDF-56 | 217193 | J033247.60$-$274452.0 | 6 | …
ALMA-HUDF-57 | … | … | 20 | A2GS33
ALMA-HUDF-58 | … | … | 17 | A2GS38
ALMA-HUDF-59 | … | J033230.35$-$274501.1 | … | …
ALMA-HUDF-60 | 209962 | … | 8 | A2GS19
ALMA-HUDF-61 | 199505 | J033244.59$-$274835.8 | 19 | …
ALMA-HUDF-62 | … | J033229.23$-$274509.7 | 11 | A2GS23
ALMA-HUDF-63 | 201793 | … | … | …
ALMA-HUDF-64 | … | … | … | …
ALMA-HUDF-65 | 214839 | … | 9 | A2GS40
ALMA-HUDF-66 | 195280 | J033240.74$-$274926.1 | … | …
ALMA-HUDF-67 | 199996 | J033228.82$-$274829.7 | 44 | …
ALMA-HUDF-68 | 198545 | J033238.73$-$274839.8 | 29 | …
ALMA-HUDF-69 | 198451 | J033243.67$-$274850.8 | 26 | …
ALMA-HUDF-70 | … | … | … | …
ALMA-HUDF-71 | … | … | … | …
ALMA-HUDF-72 | 200418 | J033241.75$-$274824.9 | … | …

##  Appendix B Source cutouts

Here we provide 5 arcsec$\,{\times}\,$5 arcsec cutouts of our ALMA sources
overlaid on JWST F356W images. In each panel of Fig. 16 we show ALMA contours
from our deep map in steps of 2${}^{n}\sigma$, with $n\,{=}\,1$, 2, 3, 4, 5
and 6. ALMA centroids are indicated by magenta circles and the positions of
JWST counterparts are indicated by green crosses. The JWST images were lightly
smoothed by a Gaussian filter with a FWHM of 1.8 pixels for presentation
purposes. Fig. 17 shows the same thing but for our wider shallow map; note
that some of the ALMA sources found in this wider map are outside of the JADES
footprint, so we leave their backgrounds blank.

Figure 16: 5 arcsec$\,{\times}\,$5 arcsec cutouts. The background shows JWST
F356W imaging, while the contours are from our deep ALMA map in steps of
2${}^{n}\sigma$, with $n\,{=}\,1$, 2, 3, 4, 5 and 6. ALMA centroids are
indicated by magenta circles and the positions of JWST counterparts are
indicated by green crosses. ALMA IDs 28, 39, and 44 all have a F356W flux
density signal-to-noise ratios less than 5, and we do not use their photometry
to fit SEDs. Figure 17: Same as Fig. 16, but for our sources in the wider map.
For ALMA sources found outside the JADES footprint the background is blank.