A Census of the Deep Radio Sky with the VLA I:
10 GHz Survey of the GOODS-N field

Abstract

We present the first high-resolution, high-frequency radio continuum survey that fully maps an extragalactic deep field: the 10 GHz survey of the Great Observatories Origins Deep Survey-North (GOODS-N) field. This is a Large Program of the Karl G. Jansky Very Large Array that allocated 380 hours of observations using the X-band ($8-12$ GHz) receivers, leading to a 10 GHz mosaic of the GOODS-field with an average rms noise $\sigma_{\rm n}=671\,\rm nJy\,beam^{-1}$ and angular resolution $\theta_{1/2}=0\hbox{$.\!\!^{\prime\prime}$}22$ across 297 $\rm arcmin^{2}$. To maximize the brightness sensitivity we also produce a low-resolution mosaic with $\theta_{1/2}=1\hbox{$.\!\!^{\prime\prime}$}0$ and $\sigma_{\rm n}=968\,\rm nJy\,beam^{-1}$, from which we derive our master catalog containing 256 radio sources detected with peak signal-to-noise ratio $\geq 5$. Radio source size and flux density estimates from the high-resolution mosaic are provided in the master catalog as well. The total fraction of spurious sources in the catalog is 0.75%. Monte Carlo simulations are performed to derive completeness corrections of the catalog. We find that the 10 GHz radio source counts in the GOODS-N field agree, in general, with predictions from numerical simulations/models and expectations from 1.4 and 3 GHz radio counts.

1 Introduction

Extragalactic surveys are essential tools to carry out statistical analysis of galaxy populations and investigate the physical processes regulating galaxy evolution throughout cosmic time. Radio continuum surveys at $\sim$GHz frequencies are particularly valuable, because they allow us to probe non-thermal processes like synchrotron emission from supernova remnants (e.g., Dubner & Giacani, 2015) and relativistic jets powered by accreting supermassive black holes (e.g., Miley, 1980). Further, thermal (free-free) radiation from Hii regions is detectable in the radio regime and dominates the total radio emission of star-forming galaxies (SFGs) at frequencies $30\,\rm GHz\lesssim\nu\lesssim 100\,\rm GHz$ (e.g., Murphy et al., 2018a; Klein et al., 2018). Radio continuum surveys, therefore, provide a unique window into the SFGs and Active Galactic Nuclei (AGN) populations (see Condon, 1992; Tadhunter, 2016, for a review). This has motivated the implementation of increasingly wider and deeper extragalactic radio surveys during the past decades (see left panel of Figure 1, and Simpson, 2017, for a review). Because the primary beam areas $(\Omega_{\rm beam}$) in radio observations are inversely proportional to the observed frequency $(\Omega_{\rm beam}\propto\nu^{-2}$), most surveys of the extragalactic radio sky have been obtained at $\approx 1-3$ GHz (e.g., Afonso et al., 2001; Seymour et al., 2004; Schinnerer et al., 2007; Ibar et al., 2009; White et al., 2012; Smolčić et al., 2017a; Heywood et al., 2021; Best et al., 2023; Hale et al., 2021, 2023). Moreover, galaxies are easier to detect at low frequencies (i.e., $\sim$GHz) where synchrotron-dominated emission leads to a steep spectral index $\alpha\sim-0.7$, which generally describes the radio spectral energy distribution (SED) of SFGs and AGN following $S\propto\nu^{\alpha}$ (e.g., Tabatabaei et al., 2017; Tisanić, K. et al., 2020; An et al., 2024).

Figure 1: Left: A summary of past, present, and future extragalactic radio surveys in the sensitivity vs survey area plane. Figure adapted from Smolčić et al. (2017a). The data are color-coded according to the observed frequency. The Square Kilometer Array Phase 1 (SKA1) data points are taken from the SKA1 Science Requirements. These predictions are made for a 10,000 hours survey across 1,000 degree² at 1.4 GHz and a 1,000 hours survey across 0.5 degree² at 9 GHz. The next generation Very Large Array (ngVLA) data points are derived using the ngVLA exposure calculator tool. We adopt the central frequencies of the Band 1 (Band 2) of the ngVLA (see Murphy, 2022, for a review), and assume 1,0000 (1,000) hours of telescope time to map a 1,000 (0.5) degree² region. For both telescopes, spatial resolutions of $0\hbox{$.\!\!^{\prime\prime}$}5$ and $0\hbox{$.\!\!^{\prime\prime}$}05$ are assumed when observing at 1.4/2.4 GHz and 8/9 GHz, respectively. Right: A compilation of deep extragalactic radio maps obtained with sub-arcsec resolution (Muxlow et al., 2005; Biggs & Ivison, 2008; Miettinen et al., 2015; Lindroos et al., 2016; Rujopakarn et al., 2016; Guidetti et al., 2017; Bondi et al., 2018; Cotton et al., 2018; Muxlow et al., 2020; Jiménez-Andrade et al., 2019, 2021), color-coded according to the RMS noise, which have been preferentially obtained at $1-3$ GHz. The 10 GHz survey of the GOODS-N field stands as one of the deepest ever obtained. It is also the first observational campaign that fully maps an entire extragalactic field at high frequencies, providing some of the highest angular resolutions ever achieved in deep radio maps.

Enabled by the improved broadband and wide-field imaging capabilities of modern radio interferometers like the Karl G. Jansky Very Large Array (VLA) and MeerKAT, it is now possible to explore the $\mu\rm Jy$ radio source population at $\sim$GHz frequencies across $\gtrsim\rm 1\,deg^{2}$ regions (e.g., Smolčić et al., 2017a; Matthews et al., 2021a; Hale et al., 2023) and, thereby, carry out systematic studies of radio-selected SFGs and AGNs out to $z\approx 5$ (e.g., Smolčić et al., 2017b; Delvecchio et al., 2017; Novak et al., 2017; Vardoulaki et al., 2019; Leslie et al., 2020; Matthews et al., 2021b; Amarantidis et al., 2023). Due to the large areal coverage of current $\sim$GHz deep radio surveys, constraints on the sub-mJy radio source counts are less influenced by sample/cosmic variance, which is essential for testing and refining theoretical models of galaxy evolution (e.g., Mancuso et al., 2015, 2017; Bonaldi et al., 2019, and references therein). Furthermore, a small fraction of $\sim$GHz radio surveys have even reached sub-arcsec angular resolutions (right panel of Figure 1; e.g., Smolčić et al., 2017a; Muxlow et al., 2020), allowing us to investigate the radio morphological properties of high-redshift compact sources (e.g., Bondi et al., 2018; Cotton et al., 2018; Jiménez-Andrade et al., 2019, 2021; Vardoulaki et al., 2019).

Despite the rapidly growing number of extragalactic radio continuum surveys, the high-frequency ($\sim 10-100\,\rm GHz$) radio sky has been sparsely explored (right panel of Figure 1). High-frequency radio surveys are important to investigate, for example:

• mechanisms for cosmic-ray energy losses,

• young radio sources whose radio spectra peak at progressively higher frequencies,

• corrections for astrophysical foregrounds in Cosmic Microwave Background (CMB) maps (e.g., de Zotti et al., 2005),

• and anomalous microwave emission (AME; e.g., Murphy et al., 2018b) arising from spinning and magnetized ultra-small dust grains.

Most important for studies on star formation, and the scope of this manuscript, high-frequency observations are sensitive to free-free emission that is a better dust-unbiased tracer of “current” star formation (e.g., Murphy et al., 2011), as opposed to synchrotron that traces cumulative history of star formation. These science topics have motivated high-frequency extragalactic surveys at $\approx 10-20$ GHz (Bolton et al., 2004; Sadler et al., 2006; Whittam et al., 2016; Huynh et al., 2019) and even 95 GHz (Sadler et al., 2008; González-López et al., 2019); nevertheless, most of the these surveys reached depths $\gtrsim 0.1\,\rm mJy$ where the dominant radio source population are AGNs. In a pioneering effort to probe the high-frequency extragalactic sky at $\mu$Jy levels, Richards et al. (1998) and Fomalont et al. (2002) carried out single-pointing VLA observations at $\approx 8.5$ GHz that reached up to a $\approx 1.5\,\rm\mu Jy\,beam^{-1}$ sensitivity and synthesized beam with FWHM of $3\hbox{$.\!\!^{\prime\prime}$}5$, which helped demonstrating that there is an increasing contribution from SFGs to the total radio source population in the $\mu$Jy regime. More recently, Algera et al. (2021) and van der Vlugt et al. (2021) obtained single-pointing VLA continuum observations at 34 GHz and 10 GHz down to a rms noise of $1.3\,\rm\mu Jy\,beam^{-1}$ and $0.41\,\rm\mu Jy\,beam^{-1}$, respectively, and angular resolutions $\gtrsim 2\hbox{$.\!\!^{\prime\prime}$}0$. These ultra-deep observations allowed Algera et al. (2021) to verify the robustness of free–free emission as a SFR indicator at high redshift; consequently, high-frequency radio emission has been used to derive the first constraints on the cosmic star formation history from free–free radio emission, which agrees with the ones inferred from other widely tested SFR indicators (Algera et al., 2022). Moreover, the aforementioned 10 and 34 GHz VLA observations led to some of the first constraints on the radio source counts in the $\mu\rm Jy$ regime at high-frequencies (Algera et al., 2021; van der Vlugt et al., 2021), albeit such results are potentially affected by sample and/or cosmic variance due to the small areal coverage of these single-pointing VLA observations.

1.1 A VLA 10 GHz Large Program in GOODS-N

While deep, high-frequency radio observations are becoming increasingly available, to date, these are limited to single-pointing maps with coarse angular resolutions. To demonstrate the feasibility of a large survey of the high-frequency radio sky with the VLA, Murphy et al. (2017) carried out a pilot program using a single-pointing in the Great Observatories Origins Deep Survey-North (GOODS-N; Dickinson et al., 2003; Giavalisco et al., 2004) at 10 GHz. Observing at this frequency has the distinct advantage of yielding sub-arcsec angular resolution imaging while probing higher rest-frame frequencies of galaxies with increasing redshift, where emission becomes dominated by thermal (free-free) radiation and directly provides a dust-unbiased measurement of massive star formation activity. By targeting the GOODS-N field one also maximizes the impact of galaxy formation and evolution research, as this is one of the best-studied extragalactic fields at optical/near-infrared wavelengths. Ancillary data available in this field include extremely deep observations from the Hubble Space Telescope (HST), James Webb Space Telescope (JWST), Spitzer Space Telescope, Chandra X-ray Observatory, Herschel Space Observatory, and the XMM-Newton Observatory (see Barro et al., 2019; Eisenstein et al., 2023; Oesch et al., 2023, and references therein). There is deep, high-resolution radio imaging in this field at 1.5 GHz (Morrison et al., 2010; Owen, 2018; Muxlow et al., 2020) and 3 GHz (Jiménez-Andrade, et al. in prep). Additional, yet shallower, 5, 5.5, and 8 GHz data are also available for a small fraction of the GOODS-N field (Richards et al., 1999; Guidetti et al., 2017; Gim et al., 2019).

The pilot program of Murphy et al. (2017) proved that combining multi-configuration VLA 10 GHz data significantly improves the capability to recover integrated flux densities of both extended and compact sources, measure source sizes, and obtain radio spectral indices and thermal fractions using the existing radio imaging in GOODS-N. Combining information from 10 GHz images with circular synthesized beams with FWHM of $1\hbox{$.\!\!^{\prime\prime}$}0$ and $0\hbox{$.\!\!^{\prime\prime}$}22$ and rms noises of $1.1\rm\,\mu Jy\,beam^{-1}$ and $572\rm\,nJy\,beam^{-1}$, respectively, Murphy et al. (2017) report the detection of 38 radio sources (above the 3.5$\sigma$ level) with an optical and/or near-infrared counterpart with a median redshift of $1.24\pm 0.25$. The resolution of $0\hbox{$.\!\!^{\prime\prime}$}22$ sufficed to derive the deconvolved FHWM of all the 32 radio sources detected in the high-resolution map, leading to a median effective radius of $69\pm 13$ mas that translates into $\approx 509\pm 114$ pc at the median redshift of this galaxy sample. These radio sizes are a factor $\sim 7$ smaller, on average, than the optical size, suggesting that star formation is centrally concentrated in these 10 GHz-detected galaxies at $z\approx 1.24$.

Motivated by the results from the pilot program reported in Murphy et al. (2017), we have conducted a VLA Large Program to produce a deep, high-resolution mosaic of the entire GOODS-N field at 10 GHz. The new set of observations consists of 17 VLA pointings with an angular resolution and sensitivity similar to that obtained in our single-pointing pilot program (Murphy et al., 2017). As a result, this is the first observational campaign that fully maps an entire extragalactic field at high frequencies and high angular resolutions reaching sub-$\mu$Jy sensitivities. Specifically, a deep 10 GHz mosaic with an angular resolution of $0\hbox{$.\!\!^{\prime\prime}$}22$ is necessary to probe the spatial distribution of massive, dust-obscured star formation in galaxies at $0.5\lesssim z\lesssim 4$. Measuring the structure in the radio regime relative to the optical/ultraviolet, for example, will be key to link the level and nature of star formation and AGN activity to the stellar mass buildup in galaxies.

Here, we report the radio continuum data products (mosaics and radio source catalogs) and inferred 10 GHz radio source counts. This is the first of a series of manuscripts that will explore the radio source populations in the GOODS-N field using the 10 GHz data reported here and recently obtained deep, high-resolution 3 GHz data (Jiménez-Andrade, et al. in prep.).

This manuscript is organized as follows. In Section 2, we describe the 10 GHz VLA data set and the imaging procedure. Section 3 reports the source extraction and properties of the radio source catalogs. We assess the reliability of the radio source catalog in Section 4, while the inferred 10 GHz radio source counts are presented in Section 5. A summary and conclusions from this work are given in Section 6.

2 Observations, data reduction, and imaging

2.1 Very Large Array Observations

A total of 380 hours of observations were taken from September 2016 to March 2018 with the VLA towards the GOODS-N field using the X-band ($8-12$ GHz) receivers (Project code: VLA 16B-320; Principal Investigator: Eric J. Murphy). The data cover a bandwidth of 4096 MHz, separated into 32 128 MHz-wide spectral windows (SPWs), and are centered at 10 GHz. The observations were obtained with a 3s signal-averaging time and full polarisation mode, albeit this manuscript only reports the total intensity mosaic and associated radio source catalog.

300 hours of observations were taken in the A-configuration (with a maximum baseline $B_{\rm max}=36.4$ km) of the VLA to provide the best possible angular resolution to resolve the radio emission of high-redshift SFGs. These data are complemented with 65 and 15 hours of observations in the B and C configuration (with $B_{\rm max}=11.1$ and 3.4 km), respectively, to improve our sensitivity to low surface brightness structures.

Figure 2: Pointing pattern used to obtain the 10 GHz mosaic of the GOODS-N field. The pointing centers are shown by the crosses. The circles extend out to a radius of $3\farcm 82$ where the VLA primary beam response is 0.1. The blue outline illustrates the footprint of the HST/WFC3 F160W imaging from the CANDELS program (Grogin et al., 2011; Koekemoer et al., 2011).

To obtain a nearly uniform sensitivity across the GOODS-N field, seventeen pointings were chosen to obtain a hexagonal-pattern mosaic (see Figure 2). The separation between the pointings center is $\rm HPBW/\sqrt{2}=3\farcm 18$, where the half-power beam width (HPBW) of the VLA at 10 GHz is $\approx 4\farcm 5$. The seventeen pointings were observed during each of the 5 hours-long observing runs used during the observations. At the beginning of each run, 3C 286 was observed during $\approx$15 mins for flux density scale, polarization angle, and bandpass calibration. Then, J1302+5748 was observed for gain and phase calibration during $\approx 1.5$ min every $\approx 15$ min when using the A and B configuration, or every $\approx 30$ min when observing with the C configuration. Each pointing was visited once during the observing run and observed for $\approx 15-20$ min.

2.2 Data Calibration

We used the VLA calibration pipeline (version 5.6.2-3), implemented in the Common Astronomy Software Applications (CASA; McMullin et al., 2007; CASA Team et al., 2022) package, to process the 76 scheduling blocks from our data set and obtain calibrated measurement sets (MSs). This pipeline is optimized to work for Stokes I continuum data by performing basic flagging (e.g., shadowed data, edge channels of sub-bands, radio frequency interference) and deriving/applying delay, bandpass, and gain/phase calibrations. We used the pipeline “weblog” to verify the quality assurance (QA) of each flagging and calibration step for all the calibrated MSs. Additionally, to further evaluate the pipeline results, we imaged the 17 pointings per each MS and inspected the resulting 1292 images. This QA process led to the identification of defective scans arising from bad weather conditions (i.e., high phase rms values from the Atmospheric Phase Interferometer), low elevations, and high wind speeds during the observations taken in the A configuration. These defective scans, that correspond to only 3.2% of the total data, were discarded outright. We also inspected the “amplitude vs frequency” diagnostic plots in the “weblog” and found a satisfactory performance of the pipeline in flagging radio frequency interference (RFI). We verified that additional flagging of RFI remaining from the pipeline has a negligible impact on the imaging quality. Finally, we split the 76 calibrated MSs into 17 MSs containing all the data from the pointings/fields (i.e., A, B, and C configuration observations) used to cover the full GOODS-N field.

Figure 3: RMS noise across the 10 GHz mosaic of GOODS-N at $0\hbox{$.\!\!^{\prime\prime}$}22$ (left) and $1\hbox{$.\!\!^{\prime\prime}$}0$ (right) resolution. The blue outline illustrates the footprint of the HST/WFC3 F160W imaging from the CANDELS program (Grogin et al., 2011; Koekemoer et al., 2011). Due to the pointing layout, the rms noise variations among the pointing centers are less than 10%.

2.3 Imaging

The calibrated MSs with the A, B, and C configuration data from the 17 fields chosen to cover the GOODS-N field were imaged with tclean in CASA. We adopted the Multi-Term Multi Frequency Synthesis (MTMFS) imaging mode (Rau & Cornwell, 2011) that performs multi-scale and multi-term cleaning for wideband imaging. We set the number of Taylor coefficients used in the spectral model to nterms=2, i.e., the spectrum is considered as a straight line with a slope, to take into account variations of the spectral structure across the image. To reconstruct the emission of complex, extended radio sources through the multi-scale cleaning implemented in MTMFS, we look for scales extending up to 16 times the FWHM of the synthesized beam. In addition, we implement the W-projection algorithm that corrects for a non-zero w-term arising from the sky curvature and non-coplanar baselines in wide-field imaging. In practice, this algorithm hinders the presence of artifacts around sources away from the phase center. After extensive testing with several values for the number of W-projection planes to use, we find that wprojplanes=64 leads to adequate imaging quality out to the regions where the primary beam response drops to $10\%$. Self-calibration was not implemented due to the faint nature of the radio sources.

2.3.1 High-Resolution Mosaic

Following the tests performed for the pilot survey (Murphy et al., 2017), we adopt the Briggs weighting with robust=0.5. The native, synthesized beam of the combined A, B, and C configuration observations is fairly Gaussian with major and minor FWHMs $\theta_{\rm maj}=0\hbox{$.\!\!^{\prime\prime}$}23\times\theta_{\rm min}=0\hbox{$.\!\!^{\prime\prime}$}20$. For simplicity, we specified a circular Gaussian restoring beam with $\rm FWHM=0\hbox{$.\!\!^{\prime\prime}$}22$, as in our pilot survey (Murphy et al., 2017), to image each pointing individually out to a primary beam response of $10\%$. Major cycles of tclean are run in parallel with the option parallel=True and cleaning stops once the residuals are four times the rms noise. The final 17 images have $10,000\times 10,000$ pixels with a pixel scale of $0\hbox{$.\!\!^{\prime\prime}$}05$, covering a $8\farcm 33\times 8\farcm 33$ region. Following the imaging of the 17 VLA pointings towards the GOODS-N field, we used the task widebandpbcor to perform a wideband primary-beam correction. Then, the resulting 17 images are combined in a weighted fashion with the task linearmosaic to obtain a linear mosaic, $I^{lm}(x)$, given by

$$I^{lm}(x)=\frac{\sum_{p}A(x-x_{p})I_{p}(x)}{\sum_{p}A^{2}(x-x_{p})},$$ (1)

where $A$ is the VLA primary beam at 10 GHz, $I_{p}$ is $p^{th}$ deconvolved image, and $x_{p}$ the pointing center. The resulting mosaic covers a total area of $\rm 297\,arcmin^{2}$ and is centered at J2000 right ascension (RA) $12^{\rm h}36^{\rm m}55^{\rm s}$ and declination (DEC) $+62^{\circ}14^{\prime}15^{\prime\prime}$.

The distribution of the rms noise across the mosaic is shown in Figure 3. We reach a point source sensitivity of $645\,\rm nJy\,beam^{-1}$ at the pointing centers, with noise variations among these centers less than 10%. As also observed in the cumulative distribution of area vs rms noise level (Figure 4), the sensitivity remains nearly homogeneous within the central $\approx\,120\,\rm arcmin^{2}$ region. Our mosaic extends beyond the area covered by the HST imaging of $170\,\rm arcmin^{2}$, albeit the sensitivity in such outer regions ranges from $\rm 1-6\,\mu Jy\,beam^{-1}$.

Finally, we produced a non-primary-beam corrected (“flat noise”) mosaic by reverting the weights (on a pixel-by-pixel basis) used to generate the primary-beam-corrected mosaic with Equation 1. This mosaic facilitates the source extraction procedure and Monte Carlo simulations (see Section 3 and 4) as it prevents the presence of noisy edges. Likewise, it allows us to inspect the pixel brightness distribution in the mosaic without being affected by the primary beam attenuation. As observed in the left panel of Figure 5, the noise amplitude distribution is fairly Gaussian with a clear excess of pixels with flux density above five times the rms noise $\sigma_{\rm n}=\rm\,671\,nJy\,beam^{-1}$.

Figure 4: Cumulative distribution of the total area covered down to a given rms noise level. The black and gray lines correspond to the distribution of the rms noise level at $0\hbox{$.\!\!^{\prime\prime}$}22$ and $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution, respectively, within the footprint of the HST/WFC3 F160W imaging (covering 170 arcmin²). The orange and red lines show the distribution of the rms noise at $0\hbox{$.\!\!^{\prime\prime}$}22$ and $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution, respectively, across the entire 10 GHz mosaic of GOODS-N (covering 297 arcmin²). Since our VLA data extend beyond the area covered by existing HST images, the rms noise remains nearly homogeneous within the central $\approx$120 arcmin².

In the following, we refer to the mosaic with a beam FWHM $0\hbox{$.\!\!^{\prime\prime}$}22$ as our high-resolution mosaic. This mosaic is fundamental to derive the structural measurements of high-redshift radio sources in the GOODS-N field.

2.3.2 Low-Resolution Mosaic

We produced a low-resolution, ($u,v$)-tapered mosaic with a $1\hbox{$.\!\!^{\prime\prime}$}0$ synthesized beam following the same approach as in Section 2.3.1, i.e., we imaged each pointing individually with tclean using a pixel scale of $0.25\,\rm arcsec/pixel$, applied primary beam corrections with widebandpbcor, and combined the deconvolved images with linearmosaic to get the low-resolution mosaic covering the same 297$\rm\,arcmin^{2}$ region as in the high-resolution one. In this case, after extensive testing, we adopt robust=2.0 to minimize the noise level. The point source sensitivity at the pointing centers of the $1\hbox{$.\!\!^{\prime\prime}$}0$ tapered mosaic is $\approx 920\,\rm nJy\,beam^{-1}$ (right panel of Figure 3) and, similar to the high-resolution mosaic, the rms noise across the low-resolution mosaic is nearly homogeneous within the central 120$\rm\,arcmin^{2}$ (see Figure 3 and 4).

We also obtained a non-primary beam corrected (“flat noise”) version of this low-resolution mosaic (see Figure 6). Its pixel brightness distribution is accurately described by a Gaussian function with $\sigma_{\rm n}=0.968\,\rm\mu Jy\,beam^{-1}$ (see right panel of Figure 5). Besides, the distribution exhibits a clear excess of pixels with flux density values above five (and even three) times the noise amplitude, as well as the absence of pixels with negative values beyond the expected Gaussian distribution.

While the rms noise of the low-resolution mosaic is 44% higher than the high-resolution one, the corresponding brightness temperature rms of the low-resolution mosaic is 19 times lower than its high-resolution counterpart. As a result, the $1\hbox{$.\!\!^{\prime\prime}$}0$ tapered mosaic allows us to better probe the faint and extended radio emission of high-redshift radio sources. Moreover, as detailed in Section 3, the low-resolution map leads to a significantly lower fraction of spurious sources, rendering this $1\hbox{$.\!\!^{\prime\prime}$}0$ tapered mosaic the preferred one to perform our blind source extraction and obtain the master radio source catalog.

Figure 5: Pixel brightness distribution of the 10 GHz mosaic of the GOODS-N field at $0\hbox{$.\!\!^{\prime\prime}$}22$ (left) and $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution (right). The pixel values shown here are uncorrected for primary-beam attenuation. The gray regions are Gaussian fits to the observed distributions. The rms noise values are shown in the upper-right corners, while the $+/-5\sigma_{\rm n}$ levels are illustrated by the dashed vertical lines. The noise of both versions of the 10 GHz mosaic of the GOODS-N field follows a Gaussian probability distribution, indicating that the mosaics are not limited by dynamic range.

Figure 6: The 10 GHz mosaic of the GOODS-N field covering 297 arcmin². Here, we show the non-primary-beam-corrected version of the mosaic at $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution with a median rms noise of 0.968 $\rm\mu Jy\,beam^{-1}$. The red circles highlight the position of the 256 radio sources in our master 10 GHz catalog of the GOODS-N field.

2.3.3 A Note About Joint Deconvolution and the A-Term Correction

Given the large data volume of the survey, totaling $\approx$ 40 TB of calibrated MSs, joint deconvolution was unfeasible even with parallel processing recently implemented in mpicasa and the computing resources available to us at NRAO. Likewise, we attempted correcting for the A-term to take into account the baseline, time, and frequency dependence of the aperture illumination pattern (AIP) of the antennas (Bhatnagar et al., 2013). This can be done with the gridder “awproject” in tclean that applies the A- and W-Projection algorithms. Nevertheless, the computing cost of the AW-Projection is significantly larger than standard imaging –even with parallelization. Hence, joint deconvolution of our entire data set with the AW-projection was nonviable.

To verify that the adopted approach to image our data set (i.e., imaging each pointing individually and combining them with linearmosaic) leads to a mosaic that is consistent with the ones produced via joint deconvolution and the AW projection, we perform the following tests. We downscale our data set and image only four out of the 17 pointings in our survey via joint deconvolution with the gridder “mosaic” and “awprojection”. We find no significant difference between the pixel brightness distribution, rms noise, number of detected sources, fraction of spurious sources, and presence of imaging artifacts between the maps produced via joint deconvolution with gridder “mosaic”$/$“awprojection” and the map obtained with our adopted approach. Moreover, we find that the structural parameters of detected sources (integrated flux density and major FWHM) in the three maps differ, in general, by $\lesssim 10\%$.

Considering that the next generation VLA will be regularly producing radio surveys with ten times better angular resolutions and sensitivity than the 10 GHz survey of GOODS-N (Murphy, 2022), it is worth stressing that massive computing resources will be needed to process such amount of data. Ongoing efforts to analyze and optimize the size-of-computing for ngVLA synthesis imaging are being taken, concluding that parallelization and implementations based on Graphics Processing Units (GPUs) and Field Programmable Gate Arrays (FPGAs) have the potential to reduce the computing costs of the next generation VLA. For example, an experimental project led by the NRAO has shown that a nationwide grid of computers with GPUs can reduce the imaging running time by two orders of magnitude (S. Bhatnagar priv. comm.).

3 Catalog

We adopt the low-resolution mosaic (with robust=2 and uv-tapered to a $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution) to derive our master radio source catalog based on the following considerations.

First, after extensive tests, we find that a uv-tapered map with a $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution leads to the lowest fraction of spurious sources compared to (Section 4.3) any other map with higher resolution. Note that the pixel brightness distribution of the high-resolution mosaic (left panel of Figure 5) shows a tail of negative sources that deviates from the Gaussian model. These spurious sources are randomly distributed across the high-resolution mosaic and disappear in the low-resolution mosaic (right panel of Figure 5). Moreover, since using robust=2 leads to the best sensitivity, the low-resolution map is the best alternative to increase the number of detections while minimizing the presence of spurious sources. A detailed description of the tests performed to unveil the dependence between the number of spurious sources and angular resolution of the map will be part of an upcoming technical report/manuscript.

Second, a master catalog from the low-resolution mosaic simplifies the implementation of Monte Carlo simulations to derive completeness and flux density boosting corrections (Section 4.1 and 4.2). Specifically, mock radio sources “observed” at $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution can be generated with a single 2D Gaussian model, instead of more complex models needed to reproduce the compounded radio structures that are revealed in the high-resolution mosaic with a $0\hbox{$.\!\!^{\prime\prime}$}22$ resolution (see Figure 8).

3.1 Source Extraction

We use the Python Blob Detector and Source Finder (PyBDSF; Mohan & Rafferty, 2015) to obtain radio source catalogs. The source extraction is performed in the “flat noise” mosaic to mitigate the effect of the noise edges on the calculation of the rms noise map, which is derived using the PyBDSF´s suggested values for the box size and step size. We adopt a threshold to identify the islands of contiguous emission (thresh_isl) of 3$\sigma$, where $\sigma$ is the local rms noise. The source detection threshold (thresh_pix) is set to 5$\sigma$. After visually inspecting the resulting 266 catalog entries, we find and remove three entries linked to artifacts around a bright radio source with $S_{\rm P}\approx\rm 410\,\mu Jy\,\rm beam^{-1}$. Also, we group 13 catalog entries into 6 multi-component, extended, and complex sources (see Section 3.2). Our master catalog, therefore, comprises 256 radio sources detected in the low-resolution mosaic (see Figure 7).

Figure 7: A flowchart that outlines the process to obtain our master catalog from the low-resolution mosaic. The flowchart also illustrates the fraction of sources that have counterparts in the high-resolution mosaic and the fraction of multi-component sources in our master catalog.

In addition, to provide more robust information on the major and minor FWHM of these radio sources, we run PyBDSF on the “flat noise” high-resolution mosaic using the detection thresholds and rms map derivation procedure adopted for our master catalog. By matching the catalogs from high- and low-resolution mosaics (using a $1\hbox{$.\!\!^{\prime\prime}$}0$ radius) and visually inspecting the radio sources, we find the following (see Figure 7).

179 sources from our master catalog are detected with a peak signal-to-noise ratio SNR$\geq$5 in both the low- and high-resolution mosaic, of which 168 appear as compact sources in both mosaics, and 11 appear as multi-component radio sources in one or both of the mosaics.

77 sources from our master catalog are only detected in the low-resolution mosaic, of which 76 are compact sources, and one is a multi-component source. In Figure 8, we present examples of the four types of radio sources in our master catalog listed above.

Lastly, in Figure 9 we focus on the 168 compact radio sources detected in the low- and high-resolution mosaics and compare their integrated flux densities measured at both angular resolutions. In general, the flux densities derived at different resolutions are similar. At the faint end ($S_{\rm I}\lesssim\rm 10\,\mu Jy$), however, the low-resolution imaging allows us to retrieve a factor 1.5 more emission –on average– than that observed in the high-resolution mosaic.

3.2 Measuring Flux Densities of Multi-Component Radio Sources

To improve the flux density measurements of the multi-component radio sources in our mosaics, we run PyBDSF in the interactive mode and follow the software´s recommendations to fit extended and complex radio sources. To improve the island determination, we set rms_map=False and mean_map=‘const’ to use a constant mean and rms value across the fits file cutouts containing the multi-component radio sources. We also set flag_maxsize_bm=50 to fit larger Gaussian components when necessary, and atrous_do=True to fit Gaussians to the residual image and model any extended emission missed in the standard fitting. Finally, we adjust the threshpix parameter to improve the fitting and, if possible, to force PyBDSF to associate multiple Gaussian components into a single, extended radio source.

Figure 8: Stamps from the 10 GHz mosaic of the GOODS-N field at $1\hbox{$.\!\!^{\prime\prime}$}0$ (left panels) and $0\hbox{$.\!\!^{\prime\prime}$}22$ arcsec resolution (right panels). The blue rectangle in the low-resolution stamp highlights the zoomed-in region displayed in the stamp at $0\hbox{$.\!\!^{\prime\prime}$}22$ resolution. These examples illustrate the four types of radio sources in our master catalog from the low-resolution mosaic (see Figure 7): (first/top row) 168 compact radio sources with a single counterpart detected in the high-resolution mosaic, (second row) 11 sources that appear as multi-component in the low- and/or high-resolution mosaic, (third row) 76 compact radio sources without a counterpart detected in the high-resolution mosaic, (fourth/bottom row) 1 multi-component source without a counterpart detected in the high-resolution mosaic. Contour levels are at 3, 5, and 8 times the rms noise.

3.3 Radio Size Estimates

A total of 168 compact sources from our master catalog have a single counterpart in the high-resolution mosaic. Radio size estimates for these 168 sources are derived from the high-resolution mosaic. In Table A, we report the deconvolved FWHM ($\theta$) along the major and minor axis of the radio sources which, in the case of a circular beam, are given by:

$$\theta=\left(\phi^{2}-\theta_{1/2}^{2}\right)^{1/2},$$ (2)

where $\phi$ is the FWHM of the fitted major or minor axis of the source and $\theta_{1/2}$ is the FWHM of the synthesized beam. The uncertainties on the deconvolved FWHM that PyBDSF reports are the same as the uncertainties on the FWHM values prior deconvolution. That seriously underestimates $\sigma_{\theta}$ for marginally resolved sources, so we estimate $\sigma_{\theta}$ using Equation (3) from Murphy et al. (2017) instead:

$$\left(\frac{\sigma_{\theta}}{\sigma_{\phi}}\right)=\left[1-\left(\frac{\theta_{1/2}}{\phi}\right)^{2}\right]^{-1/2},$$ (3)

where $\sigma_{\phi}$ is the uncertainty on the fitted FWHM prior to deconvolution. In some cases, PyBDSF reports unrealistic $\phi$ values (i.e., fitted FWHM equal to or smaller than the synthesized beam). The corresponding deconvolved FWHM are, therefore, reported as 0 in Tables A and A, while the associated uncertainty corresponds to the error on the fitted FWHM ($\sigma_{\phi}$).

Following Murphy et al. (2017), we deem sources meeting the criterion $\phi_{\rm M}-\theta_{1/2}\geq 2\sigma_{\phi_{\rm M}}$ as confidently resolved along their major (M) axis. Out of the 168 sources in our master catalog with a single counterpart in the high-resolution mosaic, 92 are confidently resolved.

In the case of sources that are not confidently resolved along their major (nor minor) axis, the peak brightness value approaches that of the integrated flux density. Therefore, we estimate the geometric mean of the peak brightness and integrated flux densities and adopt the resulting value as the best estimate for the integrated flux density ($S_{*}$; reported in Tables A and A as well). For sources whose major axes are resolved, the best estimate for the source´s integrated flux density is simply that reported by PyBDSF.

Figure 9: Comparison between the integrated flux densities of 168 radio sources detected in the high- and low-resolution mosaics of the GOODS-N field. The blue dashed line is set at 1. At $S_{\rm I}^{0.22}>\rm 10\,\mu Jy$, the flux densities at low and high resolution typically are equal within $\approx 25\%$. The red points show the median $S_{\rm I}^{1.00}/S_{\rm I}^{0.22}$ ratios within a $\approx$0.25 dex $S_{\rm I}^{0.22}$ bin. The error bars denote the 16th and 84th percentile values.

3.4 Astrometric Accuracy

We compare the positions of radio sources detected with the European VLBI Network (EVN) at 1.6 GHz (Radcliffe et al., 2018) and their counterparts in our high-resolution VLA 10 GHz mosaic of the GOODS-N field. Out of the 31 VLBI-detected sources, 26 are detected above peak $\rm SNR\gtrsim 10$ in our mosaic. The remaining five VLBI-detected sources either lie outside the footprint (four sources) or at the edge (one source) of the 10 GHz mosaic –where the sensitivity drops by a factor of ten with respect to the central region. We find a median offset between the EVN 1.6 GHz and VLA 10 GHz positions of only 2.3 milliarcsec in RA and 5.6 milliarcsec in DEC (see Figure 10). Using the positions from the low-resolution VLA 10 GHz mosaic leads to median offsets of 29.7 milliarcsec in RA and 5.7 milliarcsec in DEC. The aforementioned mean positional offsets are $\lesssim 8$ times smaller than the pixel scale of the low- and high-resolution mosaics. We did not correct the catalogs’ entries for the respective mean positional offsets derived here.

Figure 10: Comparison between the central positions of 26 radio sources in our high-resolution VLA 10 GHz mosaic and their counterparts in the EVN 1.6 GHz data from Radcliffe et al. (2018). The red star shows the median value from the offset distribution that is $\lesssim 5$ milliarcsec.

3.5 Summary of Released Catalogs

Our master source catalog comprises 256 radio sources detected in our $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution mosaic of GOODS-N, out of which 12 are multi-component. The flux densities of the 256 sources (at both angular resolutions when available) are reported in Table A. Size estimates obtained from the $0\hbox{$.\!\!^{\prime\prime}$}22$ resolution mosaic are presented as well. In Tables A and A, we report the properties of the individual components/sources in the low- and high-resolution mosaics, respectively, that are grouped into the 12 multi-component sources in our master catalog.

Henceforth, all the analyses to characterize our master catalog (Section 4) and derive the radio source counts (Section 5) are based on the 1$\hbox{$.\!\!^{\prime\prime}$}0$ resolution mosaic of GOODS-N—from which our master radio source catalog is extracted.

Figure 11: Left: Completeness of the 10 GHz source catalog of the GOODS-N field as a function of integrated flux density and SNR (for unresolved sources). The thin red lines show the completeness curve for each of the 500 mock catalogs in our Monte Carlo simulations. The solid and dashed blue lines show the 50th and 16th/84th percentiles, respectively. Right: Flux boosting as a function of SNR of mock sources in our Monte Carlo simulations. A 2D histogram is shown in the background. The darker colors indicate a higher density of sources per flux density bin. The blue markers show the 50th percentile of the input-to-output flux density ratios per 0.2 dex bin, while the error bars represent the 16th and 84th percentiles.

4 Radio source counts corrections

We assess the reliability of the master 10 GHz radio source catalog of GOODS-N by deriving corrections factors to account for the completeness, flux boosting, false detections, and resolution bias.

4.1 Completeness

To determine the number of sources that exist in a given region of the sky (above a detection limit) but are missed in our mosaic/catalog due to the adopted observational and detection procedure, we perform extensive Monte Carlo simulations as follows.

1. We infer the probability distribution function (PDF) of the peak brightness ($S_{\rm P}$) and deconvolved major/minor FWHM ($\theta_{\rm M,m}$) of the sources in our catalog by fitting an exponential and a half-norm function ($f(x)=\sqrt{2/\pi}\exp(-x^{2}/2)$; for $x\geq 0$), respectively.

2. We use the inferred PDFs of the $S_{\rm P}$ and $\theta_{\rm M,m}$ to generate a mock sample of radio sources. We verify that the distribution of the integrated flux density from the mock catalog matches, in general, that of the observed catalog.

3. 1000 mock radio sources are injected in the mosaic at random positions and position angles, under the condition that mock sources are located $\approx 9\,\rm arcsec$ away from real or other mock sources. These mock sources follow a $S_{\rm P}$ distribution that includes values as low as $3\mu\rm Jy\,beam^{-1}$, i.e., $\approx 3\times$ the average rms noise of our mosaic. By injecting these faint sources we consider the effect of the noise in boosting/decreasing the flux density of sources with $\rm SNR\approx 5$.

4. We take into account the rms noise variations in our mosaic, mainly arising from the primary beam attenuation of the 17 pointings, to derive our mock catalogs and completeness corrections as follows. Let us first gauge an illustrative case of a compact source with $S_{\rm I}\approx 50\,\rm\mu Jy$ that lies at the edge of our mapped region. While compact sources with such $S_{\rm I}$ values are robustly detected with a $\rm SNR\approx 50$ in the central region of our mosaic, a source with $S_{\rm I}\approx\rm 50\mu Jy$ is detected with $\rm SNR\approx 5$ at the outskirts of the map where the primary beam response drops to 10%, which hinders the completeness of our catalog for such a hypothetical flux density value. Therefore, to consider the effect of the primary beam response in our completeness correction, the input flux density of the mock sources is reduced depending on their position in the mock mosaic.

5. We repeat steps two to four and generate 500 mock mosaics, translating into a half million sources in our Monte Carlo simulations. Then, we obtain the corresponding 500 mock catalogs with PyBDSF by applying the same detection parameter criteria used to obtain our 10 GHz catalog of the GOODS-N field.

6. We derive the completeness of our catalog by comparing the number of detected sources with the number of injected sources, per integrated flux density bin, for all the mock mosaics/catalogs. The resulting 500 completeness curves, already corrected by the primary beam attenuation across the mosaic, are shown in the left panel of Figure 11. We adopt the median trend and the 16th/84th percentiles as our best completeness values ($C_{\rm comp}$) and associated errors, which are also reported in Table 1.

Table 1: Completeness correction factors ($C_{\rm comp}$) for the GOODS-N 10 GHz catalog as a function of integrated flux density.

Integrated flux density	$C_{\rm comp}$	Error
($\mu\rm Jy$)
$\leq$3.4	0.00	0.00
4.3	0.09	0.05
5.4	0.18	0.06
6.8	0.30	0.06
8.5	0.42	0.06
10.7	0.53	0.05
13.5	0.64	0.05
17.0	0.72	0.05
21.4	0.78	0.04
26.9	0.84	0.04
33.9	0.89	0.03
42.7	0.92	0.03
53.7	0.95	0.03
67.6	0.97	0.03
85.1	0.98	0.02
107.1	1.00	0.00
134.9	1.00	0.00
169.8	1.00	0.00
213.8	1.00	0.00
269.1	1.00	0.00
338.8	1.00	0.00
426.6	1.00	0.00

4.2 Flux Boosting

To characterize the effects of flux boosting (e.g., Coppin et al., 2005; Casey et al., 2014) in our sensitivity-limited mosaic, we contrast the input and output flux density ($S_{\rm I}^{\rm in}$ and $S_{\rm I}$, respectively) of the mock sources in our Monte Carlo simulations. We find that flux densities are boosted by $\approx 14\%$ in the lowest SNR bin centered at $S_{\rm I}\approx 6\mu\rm Jy$ (Figure 11). At SNR larger than 15, flux densities are boosted by less than 5%. Consequently, we do not apply this correction to the flux densities in the catalog.

4.3 False Detections

We determine the fraction of spurious sources in our catalog by performing the source extraction with PyBDSF on the inverted (i.e., multiplied by -1) mosaic. To this end, we use the same detection parameters adopted to obtain the 10 GHz catalog of GOODS-N. Two spurious sources are detected with $\rm SNR\approx 5$, leading to a total fraction of spurious sources in our catalog of only 0.75%. This implies a notably high “fidelity” parameter of 0.99, defined as $1-N_{\rm neg}/N_{\rm pos}$ with $N_{\rm neg}$ and $N_{\rm pos}$ the number of negative and positive detections, respectively (Decarli et al., 2020).

Comparing the number of spurious sources with the number of sources detected in our mosaic per SNR bin, we find a fraction of spurious sources of 4% within $5.0\leq\rm SNR<5.5$, and 0% for SNR larger than 5.5. These sources are detected at the outskirts of the map, where the primary beam response is 0.1414 and 0.4359, and have a total (primary beam-attenuated) integrated flux density of $6.59\pm 2.18\rm\mu Jy$ and $9.69\pm 2.90\rm\mu Jy$, respectively. Considering the primary-beam corrected flux densities, the fraction of spurious sources ($f_{ss}$) is zero in all $S_{\rm I}$ bins except that spanning from 15 to $50\rm\mu Jy$ where $f_{ss}=0.02$ (see Table 2).

4.4 Resolution Bias

Detecting sources in our SNR thresholded mosaic relies on the peak brightness. An unresolved source in our mosaic with $S_{\rm I}\approx\rm 10\mu Jy$, for example, has a greater probability of being detected than an extended source with the same $S_{\rm I}$ value but lower peak brightness. This effect will hinder the number of detections, particularly close to our detection limit. As a result, the population of extended, low surface brightness sources might be underrepresented in our original catalog and mock mosaics. We address this so-called “resolution bias” by following the analytic methodology presented in van der Vlugt et al. (2021) that is summarized in the following lines.

The resolution bias correction factor ($f_{\rm rb}$) is given by

$$f_{\rm rb}=\left[1-h(>\theta_{\rm max})\right]^{-1},$$ (4)

where $h(>\theta_{\rm max})$ is the fraction of sources expected to be larger than the maximum angular size ($\theta_{\rm max}$) that our detection procedure is sensitive to. Such a fraction can be inferred with the relation (Windhorst et al., 1990):

$$h(>\theta_{\rm max})=\exp\left[-\ln(2)\left(\frac{\theta_{\rm max}}{\theta_{\rm med}}\right)^{0.62}\right].$$ (5)

$\theta_{\rm max}$ depends on the source´s integrated flux density and is expressed as

$$\theta_{\rm max}=\left[\theta_{\rm 1/2}^{2}\times(S_{\rm I}/5\sigma)\right]^{1/2},$$ (6)

with $\theta_{\rm 1/2}$ the FWHM of our synthesized circular beam. Finally, $\theta_{\rm med}$ is the median angular size of the radio source population. Windhorst et al. (1990) propose a flux-dependent size given by $\theta_{\rm med}=2(S_{\rm 1.4\,GHz})^{0.3}$, where $S_{\rm 1.4\,GHz}$ is the flux density in millijanskys that we infer by scaling our 10 GHz measurements with a spectral index of $-0.8$. We also adopt a constant value of $\theta_{\rm med}=0\hbox{$.\!\!^{\prime\prime}$}30$ that has been specifically derived for the $\mu\rm Jy$ radio source population (Cotton et al., 2018; Bondi et al., 2018). By evaluating Equation 4 we derive $f_{\rm rb}$ for a flux-dependent and constant $\theta_{\rm med}$. Here, we adopt the mean of these two values as our best estimate for the resolution bias correction factor, which is $\approx 1.40$ for integrated flux densities $\approx 5-15\mu\rm Jy$ and less than 1.20 for $S_{\rm I}\gtrsim 100\mu\rm Jy$. A list with the $f_{\rm rb}$ values per $S_{\rm I}$ bin can be found in Table 2.

5 Radio source counts

To derive the observed Euclidean-normalized number counts, following Matthews et al. (2021a), we first estimate the quantity

$$S_{\rm I}^{2}n(S_{\rm I})=\left[\frac{1}{\Omega\ln(\Delta)}\right]\sum^{n_{\rm bin}}_{i=1}S_{\rm I}^{i},$$ (7)

where $n_{\rm bin}$ is the number of sources per $S_{\rm I}$ bin, $\Omega$ is the survey area of $297\,\rm arcmin^{2}$, and $\Delta$ the logarithmic width of the bin of $\rm dex(0.5)$. We also derive the rms statistical uncertainty in $S_{\rm I}^{2}n(S_{\rm I})$ using

$$\sigma_{\rm stat}=\left[\frac{1}{\Omega\ln(\Delta)}\right]\left(\sum^{n_{\rm bin}}_{i=1}\left(S_{\rm I}^{i}\right)^{2}\right)^{1/2},$$ (8)

which is valid for bins with $n_{\rm bin}\gg 1$. We note that the uncertainty for the brightest bin (with $n_{\rm bin}=4$) is potentially larger than the quoted/plotted value due to Poissonian fluctuations.

Then, the observed Euclidean-normalized number counts are estimated by using $S_{\rm I}^{2}n(S_{\rm I})$ multiplied by $S_{\rm I}^{1/2}$. These number counts must be corrected for completeness, fraction of spurious sources, and resolution bias. Thus, to derive the corrected Euclidean-normalized number counts we calculate

$$S_{\rm I}^{2}n(S_{\rm I})=\left[\frac{1}{\Omega\ln(\Delta)}\right]\sum^{n_{\rm bin}}_{i=1}\frac{[1-f_{\rm ss}(S_{\rm I}^{i})]\times f_{\rm rb}(S_{\rm I}^{i})}{C_{\rm comp}(S_{\rm I}^{i})},$$ (9)

where $C_{\rm comp}$, $f_{\rm ss}$, and $f_{\rm rb}$ are the aforementioned completeness, spurious sources, and resolution bias correction factors, respectively. The observed and corrected Euclidean-normalized 10 GHz radio source counts of the GOODS-N field are shown in the top panel of Figure 12. These cover the flux density range $-5.28<\log(S_{\rm I}/\rm Jy)<-2.78$ in five bins of 0.5 dex-width. In Table 2, we also report the correction factors and radio source counts per flux density bin.

We compare the 10 GHz source counts with those derived by van der Vlugt et al. (2021) using an ultra-deep, single-pointing VLA data in the COSMOS field. Their X-band VLA image reaches an rms noise level of $0.41\,\rm\mu Jy\,beam^{-1}$ at the pointing center and has an angular resolution of $2\hbox{$.\!\!^{\prime\prime}$}33\times 2\hbox{$.\!\!^{\prime\prime}$}01$. As observed in the top panel of Figure 12, the 10 GHz radio source counts derived by van der Vlugt et al. (2021) span over $-5.57<\log(S_{\rm I}/\rm Jy)<-4.10$, allowing a more direct comparison between the 10 GHz number counts derived in both studies across $-5.28<\log(S_{\rm I}/\rm Jy)<-4.10$. Within this flux density regime, the 10 GHz number counts from van der Vlugt et al. (2021) are systematically higher by a factor $\approx 1.6$ than those reported here. We discuss this discrepancy within the context of sample and cosmic variance in Section 5.2.

Figure 12: Top panel: Euclidean-normalized differential radio source counts at 10 GHz from the GOODS-N field (filled red circles), which have been corrected by completeness, fraction of spurious sources, and resolution bias. The uncorrected source counts are shown as open red circles. These are the first radio source counts obtained at high frequencies ($\nu\gtrsim 5\,\rm GHz$) from an entire extragalactic deep field. For comparison, we also present the radio source counts at 10 GHz from a single-pointing VLA map in the COSMOS field (van der Vlugt et al., 2021). To illustrate the effect of cosmic variance on single-pointing VLA imaging at 10 GHz, we present the maximum/minimum radio source counts (red dotted lines) obtained from six nonoverlapping subfields each covering 32.16 $\rm arcmin^{2}$ in our 10 GHz mosaic of GOODS-N. Middle panel: A compilation of radio source counts measured at 1.4 and 3 GHz (Smolčić et al., 2017a; Matthews et al., 2021a; van der Vlugt et al., 2021; Hale et al., 2023), including the statistical counts derived from a $P(D)$ analysis at 1.4 GHz in the DEEP2 field (shaded gray region; Matthews et al., 2021a). The 1.4 and 3 GHz radio source counts are shifted to the 10 GHz observed frame by adopting a spectral index of $-0.7$, to facilitate the comparison between 10 GHz radio source counts (van der Vlugt et al., 2021, and this work) and the expected values from the more abundant 1.4 and 3 GHz measurements. No error bars are shown here to aid the visual inspection of the plot. Bottom panel: Radio source counts measured at and/or shifted to 10 GHz compared with the simulations and model from Mancuso et al. (2017) and Bonaldi et al. (2019), who made specific predictions for the radio sky population at the $\approx$10 GHz observed frame. The dotted and dashed blue lines show the radio source counts of the AGN and SFG populations, respectively, as derived by Bonaldi et al. (2019).

5.1 Comparison with 1.4 and 3 GHz Number Counts

We compare the 10 GHz radio source counts in the GOODS-N field with the more abundant measurements obtained at 1.4 and 3 GHz in the COSMOS, XMM Large Scale Structure (XMM-LSS), and DEEP2 fields (Smolčić et al., 2017a; Matthews et al., 2021a; van der Vlugt et al., 2021; Hale et al., 2023). To this end, the 1.4 and 3 GHz radio source counts are converted to the 10 GHz observed frame assuming that the radio SED is described by a power-law: $S_{\rm\nu}/S_{\rm 10}=(\nu\,[\rm GHz]/10)^{\alpha}$, where $\alpha$ is the spectral index that here is fixed to $-0.7$. Then, the 10 GHz radio source counts ($n_{\rm 10}$) are estimated from the 1.4 and 3 GHz values ($n_{\nu}$) using $n_{\rm 10}=n_{\nu}(10/\nu\,[\rm GHz])^{1.5\alpha}$. As observed in the middle panel of Figure 12, the 10 GHz radio source counts from GOODS-N follow, in general, the trend depicted by the 10 GHz number counts inferred from the lower frequency observations. Radio source counts rise from the sub-$\mu\rm Jy$ regime, they flatten out at flux densities $10\,{\rm\mu Jy}\lesssim S_{\rm I}\lesssim 100\,\rm\mu Jy$, and then continue to rise towards the bright end. The scatter of the number counts from all these studies is, however, evident. This can be attributed to the different assumptions made to correct for the resolution bias and completeness, as well as field-to-field variations due to sample and cosmic variance (as also discussed by Smolčić et al., 2017a; van der Vlugt et al., 2021; Hale et al., 2023). Furthermore, the small discrepancies between the different trends observed in the middle panel of Figure 12 are also a consequence of the simplistic assumptions made here to scale the 1.4 and 3 GHz radio source counts to the 10 GHz observed frame, i.e., a single power-law radio SED with $\alpha=-0.7$. We note, however, that adopting flatter or steeper spectral indices alleviates the discrepancies between the scaled 1.4/3 GHz and observed 10 GHz radio source counts across different flux density regimes. For example, using $\langle\alpha^{10\,\rm GHz}_{\rm 1.4\,GHz}\rangle=-0.61$—as reported for 10 GHz radio sources in our pilot survey with and without a counterpart at 1.4 GHz (Murphy et al., 2017)—increases the normalization of the scaled 1.4 GHz source counts by $\approx 0.12$ dex. This leads to an excellent agreement between the robustly constrained 1.4 GHz source counts in DEEP2 (Matthews et al., 2021a) and our 10 GHz estimates in GOODS-N for flux densities $16\,\rm\mu Jy\lesssim S_{\rm I}\lesssim 166\,\rm\mu Jy$. However, a steeper spectral index of $\langle\alpha^{10\,\rm GHz}_{\rm 1.4\,GHz}\rangle\approx-0.8$ is needed to get a better agreement between the scaled 1.4 and 10 GHz number counts for flux densities $S_{\rm I}\lesssim 16\,\rm\mu Jy$, which is compatible with the average spectral index steeper than $-0.75$ found for 6 and 8.5 GHz-detected sources below 35$\mu$Jy (Fomalont et al., 2002; Thomson et al., 2019).

5.2 The Impact of Sample and/or Cosmic Variance on Radio Source Counts

As mentioned before, sample and/or cosmic variance could be one of many factors driving the scatter of the number counts obtained from different radio surveys. To illustrate this, we focus on a comparison between the single- and multiple-pointing 10 GHz data in the COSMOS (van der Vlugt et al., 2021) and GOODS-N field, respectively. Both data sets have been obtained with the X-band receivers of the VLA, reached similar depths and comparable angular resolutions, leaving the survey area as the main variable in our comparison.

We start by splitting the 10 GHz mosaic of GOODS-N into six nonoverlapping regions each covering $32.16\,\rm arcmin^{2}$, matching the area of the single-pointing 10 GHz map of COSMOS from van der Vlugt et al. (2021). Radio source counts from these six subfields are then obtained and corrected following Equation 9, as done for the radio source counts from the full 10 GHz mosaic of GOODS-N. The maximum and minimum number counts from these subfields, shown in the top panel of Figure 12, suggest that sample variance induces an additional scatter of $\rm\approx 0.15-0.40\,dex$ in the radio source counts measured in areas as small as $32.16\,\rm arcmin^{2}$. We note that such a scatter should be considered as a lower limit to the true cosmic variance arising from the large-scale structure of the Universe because the cosmic variance in the volume probed by the GOODS-N field is significant ($\approx 10-30\%$ at redshifts $0<z<4$, Somerville et al., 2004; Driver & Robotham, 2010). The systematically higher number counts reported in van der Vlugt et al. (2021), therefore, could be a result of sample and/or cosmic variance due to the relatively small region covered by the single-pointing VLA data. Additionally, given the two times coarser angular resolution of the COSMOS data, the resolution bias (see Section 4.4) could be also driving the higher number counts; even though we adopt the same method used by van der Vlugt et al. (2021) to correct the number counts from the potentially missing population of extended, low surface brightness sources.

5.3 Comparison with Numerical Simulations and Models

A more direct comparison between the 10 GHz radio source counts from GOODS-N can be made with the numerical simulations and models reported by Mancuso et al. (2017) and Bonaldi et al. (2019), who provide specific predictions for the radio continuum sky at $\approx$10 GHz (bottom panel of Figure 12). First, Mancuso et al. (2017) employ redshift-dependent SFR functions that are mapped into bolometric AGN luminosity functions, using deterministic evolutionary tracks for the star formation and supermassive black hole accretion in an individual galaxy, to derive differential number counts at 150 MHz, 1.4 GHz, and 10 GHz. Second, Bonaldi et al. (2019) present the Tiered Radio Extragalactic Continuum Simulation (T-RECS), which links the radio source positions to those of dark matter halos in the P-Millennium simulation (Baugh et al., 2019). Covering a 25 deg² field of view, T-RECS employs redshift-dependent AGN and SFR functions to produce a set of simulated catalogs covering the frequency range from 150 MHz to 20 GHz. The simulated 10 GHz number counts presented here are obtained by interpolating the data points from the T-RECS catalogs at 920 and 1250 MHz.

As observed in the bottom panel of Figure 12, our 10 GHz number counts at flux densities $\gtrsim 16\,\rm\mu Jy$ are consistent with Mancuso et al. (2017)’s and Bonaldi et al. (2019)’s predictions within $\approx 1.2\,\sigma$. Yet, the number counts in our faintest bin, where the dominant radio source population are SFGs (see bottom panel of Figure 12), are in better agreement with the predictions from Mancuso et al. (2017). The 10 GHz source counts at $S_{\rm I}\approx 10\rm\,\mu Jy$ from Bonaldi et al. (2019) are a factor $\approx 1.5$ higher than those reported here and the predictions from Mancuso et al. (2017), which could be a result of the assumptions made to simulate the SFG populations. While Bonaldi et al. (2019) adopts an evolving FIR-radio correlation and a modified Schechter parameterization (i.e., with two characteristic slopes) to model the SFR function, Mancuso et al. (2017) do not vary the FIR-radio ratio and employs a standard Schechter function.

Table 2: Euclidean-normalized differential 10 GHz radio source counts in the GOODS-N field.

$\Delta S_{\rm I}$	$S_{\rm I}$	$N$	$n(S_{\rm I})S_{\rm I}^{5/2}$	$f_{\rm ss}$	$f_{\rm rb}$	$C_{\rm comp}$
($\mu\rm Jy$)	($\mu\rm Jy$)		($\rm Jy^{3/2}\,sr^{-1}$)
5.24-16.56	9.31	128	$0.40^{+0.09}_{-0.13}$	0.00	$1.37^{+0.05}_{-0.03}$	$0.57^{+0.12}_{-0.18}$
16.56-52.38	29.45	92	$0.69^{+0.11}_{-0.10}$	0.02	$1.30^{+0.02}_{-0.03}$	$0.81^{+0.09}_{-0.08}$
52.38-165.63	93.14	22	$0.78\pm 0.18$	0.00	$1.21^{+0.02}_{-0.03}$	$0.98^{+0.02}_{-0.04}$
165.63-523.77	294.5	10	$1.83\pm 0.62$	0.00	$1.15^{+0.01}_{-0.02}$	1.0
523.77-1656.29	931.40	4	$4.20\pm 2.30$	0.00	$1.11^{+0.01}_{-0.02}$	1.0

.

Note. — $\Delta S_{\rm I}$ is the flux bin, centered at $S_{\rm I}$, within which the radio source counts $n(S_{\rm I})$ are estimated using $N$ number of sources per bin. The listed Euclidean-normalized differential 10 GHz radio source counts have been corrected using the $f_{\rm ss}$, $f_{\rm rb}$, and $C_{\rm comp}$ factors that account for the fraction of spurious sources, resolution bias, and completeness, respectively (see Equation 9). The quoted errors in the corrected source counts are estimated by quadratically adding the fractional errors in $f_{\rm fb}$ and $C_{\rm comp}$ and the fractional rms statistical uncertainty from Equation 8

6 Summary and conclusions

We have presented the first radio continuum survey ever obtained in an entire extragalactic deep field at high frequencies: the 10 GHz survey of GOODS-N. The main data products derived from this VLA Large Program, as well as the results from the inferred 10 GHz radio source counts, are summarized below.

• Two versions of the mosaic covering an area of $297\,\rm arcmin^{2}$ have been produced. One high-resolution mosaic with a synthesized circular beam with FWHM $0\hbox{$.\!\!^{\prime\prime}$}22$ and point-source sensitivity of $671\rm\,nJy\,beam^{-1}$, and a low-resolution, (u,v)-tappered mosaic with an angular resolution of $1\hbox{$.\!\!^{\prime\prime}$}0$ and $968\rm\,nJy\,beam^{-1}$ depth.

• We have adopted the low-resolution mosaic to obtain our master 10 GHz catalog of the GOODS-N field, which comprises 256 radio sources (detected with peak $\rm SNR\geq 5$), out of which 12 are multi-component. Size and flux density estimates from the high-resolution mosaic are reported as well.

• Monte Carlo simulations have been performed to derive the completeness of our master radio source catalog as a function of integrated flux densities. For flux densities larger than 10 (100) $\rm\mu Jy$, or $\rm SNR\approx 10$ (100) for unresolved sources, our catalog reaches a completeness of 50% (100%). The total fraction of spurious sources in our master catalog is only 0.75%.

• We have derived the 10 GHz radio source counts in the GOODS-N field. Comparing our results with the 10 GHz number counts from a single-pointing VLA image in COSMOS (van der Vlugt et al., 2021), we find that the latter are systematically higher (by a factor $\approx 1.6$). This is likely a consequence of sample and/or cosmic variance arising from the small field of view of the VLA COSMOS observations.

• The 10 GHz radio source counts in the GOODS-N field agree with the expected trend from 1.4 GHz and 3 GHz radio counts that are scaled to the 10 GHz observed frame. Nevertheless, the number counts inferred from previous studies and this work scatter across $\rm\approx 0.2-0.3\,dex$, which might be driven by the different approaches used to correct for observational biases, as well as field-to-field variations due to sample and cosmic variance.

• The 10 GHz radio source counts across ($16\lesssim S_{\rm I}/\mu{\rm Jy}\lesssim 1600$) are, in general, consistent with the predictions made by Mancuso et al. (2017) and Bonaldi et al. (2019). At the faint end ($\rm S_{\rm I}\approx 10\mu\,Jy$), the predictions from Mancuso et al. (2017) offer a better description of the 10 GHz radio source counts in GOODS-N.

Since this is the deepest and most detailed radio survey of the high-frequency radio sky, the 10 GHz mosaic of GOODS-N has the potential to address a diversity of open issues in extragalactic astronomy. The synthesized beam of $0\hbox{$.\!\!^{\prime\prime}$}22$ in the high-resolution mosaic approaches the angular resolution of HST and JWST, allowing us to explore obscured and unobscured star formation in distant galaxies at similar angular resolutions. Combining 10 GHz radio continuum imaging with available 1.4, 3, and 5 GHz data in the GOODS-N field, it will be possible to explore the radio spectro-morphological properties of $\mu\rm Jy$ sources to benchmark radio continuum emission as a robust indicator of star formation at high redshifts. These studies will be reported in forthcoming manuscripts, paving the way for the observational studies of sub-$\mu$Jy radio sources (at sub-arcsec resolutions) that will be routinely obtained with the next generation VLA in the next decade (e.g., Barger et al., 2018; Murphy, 2022; Latif et al., 2024, see Figure 1).

Appendix A 10 GHz catalogs of GOODS-N

We present a sample of the master 10 GHz catalog of GOODS-N in Table A. It reports the flux densities of sources measured in the low- and high-resolution mosaics. The radio size estimates from the high-resolution mosaic are provided as well if sources are simultaneously detected in the low- and high-resolution mosaics. Additionally, in Tables A and  A, we present the information of the multiple islands of emission that constitute the multi-component radio sources as observed in the low- and high-resolution version of the mosaic, respectively.

Table 3: A sample of the catalog of 256 radio sources in the 10 GHz mosaic of the GOODS-N field at $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution. The full catalog can be found in the online version of this manuscript and at https://science.nrao.edu/science/surveys/vla-x-gn/home.

Name	RA${}^{\tablenotemark{a}}$	DEC${}^{\tablenotemark{a}}$	$S_{\rm I}^{\rm 1.0}$ ${}^{\tablenotemark{b}}$	$S_{\rm P}^{\rm 1.0}$ ${}^{\tablenotemark{c}}$	$S_{\rm I}^{\rm 0.22}$ ${}^{\tablenotemark{b}}$	$S_{\rm P}^{\rm 0.22}$ ${}^{\tablenotemark{c}}$	$S_{*}^{\rm 0.22}$ ${}^{\tablenotemark{d}}$	$\theta_{\rm M}$	$\theta_{\rm m}$	R${}^{\tablenotemark{e}}$	T${}^{\tablenotemark{f}}$
	(degs)	(degs)	$(\rm\mu Jy)$	$(\rm\mu Jy\,beam^{-1})$	$(\rm\mu Jy)$	$(\rm\mu Jy\,beam^{-1})$	$(\rm\mu Jy)$	(mas)	(mas)
J123830.67+621821.42	189.627809 $\pm$ 0.000055	62.305950 $\pm$ 0.000031	170.80 $\pm$ 32.80	37.12 $\pm$ 5.94							C
J123821.96+621823.83	189.591517 $\pm$ 0.000028	62.306619 $\pm$ 0.000015	29.86 $\pm$ 6.09	18.15 $\pm$ 2.46	20.15 $\pm$ 2.96	18.79 $\pm$ 1.63	19.46 $\pm$ 0.67	86 $\pm$ 58	0 $\pm$ 18	0	C
J123821.78+621706.68	189.590749 $\pm$ 0.000039	62.285190 $\pm$ 0.000025	22.28 $\pm$ 6.10	12.25 $\pm$ 2.29							C
J123819.81+621839.22	189.582537 $\pm$ 0.000019	62.310894 $\pm$ 0.000025	55.38 $\pm$ 7.86	20.87 $\pm$ 2.23							C
J123819.08+621826.99	189.579516 $\pm$ 0.000014	62.307498 $\pm$ 0.000010	15.31 $\pm$ 2.91	17.24 $\pm$ 1.81	23.69 $\pm$ 3.18	15.84 $\pm$ 1.38	23.69 $\pm$ 3.18	200 $\pm$ 42	104 $\pm$ 46	1	C
J123811.83+621821.90	189.549305 $\pm$ 0.000007	62.306084 $\pm$ 0.000007	31.59 $\pm$ 2.67	24.75 $\pm$ 1.31	23.54 $\pm$ 1.57	22.35 $\pm$ 0.88	22.94 $\pm$ 0.68	0 $\pm$ 27	0 $\pm$ 8	0	C
J123810.58+621729.23	189.544078 $\pm$ 0.000042	62.291454 $\pm$ 0.000016	7.69 $\pm$ 2.37	5.88 $\pm$ 1.11	8.10 $\pm$ 1.75	5.82 $\pm$ 0.79	6.87 $\pm$ 0.69	237 $\pm$ 72	0 $\pm$ 23	0	C
J123807.76+621540.57	189.532353 $\pm$ 0.000010	62.261270 $\pm$ 0.000010	15.16 $\pm$ 2.49	17.31 $\pm$ 1.57	26.82 $\pm$ 1.68	30.95 $\pm$ 1.07	28.81 $\pm$ 0.64	0 $\pm$ 8	0 $\pm$ 7	0	C
J123807.43+621650.69	189.530948 $\pm$ 0.000039	62.280746 $\pm$ 0.000032	14.58 $\pm$ 3.77	6.16 $\pm$ 1.17							C
J123805.51+621445.55	189.522969 $\pm$ 0.000045	62.245986 $\pm$ 0.000039	31.26 $\pm$ 8.03	12.62 $\pm$ 2.38							C
J123803.66+621711.37	189.515231 $\pm$ 0.000011	62.286491 $\pm$ 0.000031	43.15 $\pm$ 4.33	11.71 $\pm$ 0.94							C
J123800.92+621336.02	189.503821 $\pm$ 0.000014	62.226672 $\pm$ 0.000020	25.10 $\pm$ 5.82	23.01 $\pm$ 3.15	42.29 $\pm$ 10.32	10.94 $\pm$ 2.17	42.29 $\pm$ 10.32	376 $\pm$ 100	369 $\pm$ 99	1	C
J123758.81+621458.32	189.495042 $\pm$ 0.000024	62.249533 $\pm$ 0.000021	19.17 $\pm$ 3.91	11.18 $\pm$ 1.54	13.99 $\pm$ 3.62	5.37 $\pm$ 1.04	13.99 $\pm$ 3.62	361 $\pm$ 106	201 $\pm$ 78	1	C
J123757.01+622059.41	189.487547 $\pm$ 0.000005	62.349835 $\pm$ 0.000005	139.62 $\pm$ 9.84	127.19 $\pm$ 5.33	136.81 $\pm$ 7.54	109.45 $\pm$ 3.74	136.81 $\pm$ 7.54	123 $\pm$ 18	97 $\pm$ 20	1	C
J123755.95+621507.54	189.483125 $\pm$ 0.000025	62.252093 $\pm$ 0.000035	19.78 $\pm$ 4.23	8.31 $\pm$ 1.30	7.45 $\pm$ 2.09	4.69 $\pm$ 0.84	5.91 $\pm$ 0.67	318 $\pm$ 110	0 $\pm$ 28	0	C
J123752.75+621628.25	189.469793 $\pm$ 0.000025	62.274515 $\pm$ 0.000023	17.07 $\pm$ 3.25	7.79 $\pm$ 1.07	15.59 $\pm$ 3.84	3.04 $\pm$ 0.63	15.59 $\pm$ 3.84	602 $\pm$ 155	318 $\pm$ 90	1	C
J123752.55+621936.96	189.468945 $\pm$ 0.000040	62.326933 $\pm$ 0.000048	23.39 $\pm$ 6.43	7.99 $\pm$ 1.69							C
J123751.23+621919.02	189.463461 $\pm$ 0.000002	62.321949 $\pm$ 0.000002	118.56 $\pm$ 2.67	104.79 $\pm$ 1.42	109.07 $\pm$ 1.86	95.90 $\pm$ 0.99	109.07 $\pm$ 1.86	93 $\pm$ 6	68 $\pm$ 8	1	C
J123750.25+621359.13	189.459382 $\pm$ 0.000022	62.233091 $\pm$ 0.000031	12.27 $\pm$ 3.34	8.21 $\pm$ 1.45	6.80 $\pm$ 1.75	6.26 $\pm$ 0.95	6.52 $\pm$ 0.67	0 $\pm$ 83	0 $\pm$ 26	0	C
J123748.15+621610.37	189.450641 $\pm$ 0.000016	62.269546 $\pm$ 0.000023	7.79 $\pm$ 1.90	6.69 $\pm$ 0.98	7.99 $\pm$ 1.64	5.32 $\pm$ 0.71	6.52 $\pm$ 0.72	201 $\pm$ 64	105 $\pm$ 70	0	C
J123747.94+621442.10	189.449761 $\pm$ 0.000019	62.245027 $\pm$ 0.000017	10.33 $\pm$ 2.17	8.40 $\pm$ 1.08	7.93 $\pm$ 1.51	6.43 $\pm$ 0.75	7.15 $\pm$ 0.70	147 $\pm$ 60	49 $\pm$ 114	0	C
J123747.08+621631.90	189.446163 $\pm$ 0.000013	62.275528 $\pm$ 0.000014	12.43 $\pm$ 2.02	10.19 $\pm$ 1.02	15.09 $\pm$ 3.00	4.35 $\pm$ 0.69	15.09 $\pm$ 3.00	457 $\pm$ 98	247 $\pm$ 63	1	C
J123746.67+621738.59	189.444478 $\pm$ 0.000001	62.294053 $\pm$ 0.000001	311.87 $\pm$ 8.21	254.71 $\pm$ 1.03	258.37 $\pm$ 1.36	230.77 $\pm$ 0.73	258.37 $\pm$ 1.36	86 $\pm$ 2	65 $\pm$ 3	1	C
J123745.88+621434.83	189.441178 $\pm$ 0.000060	62.243010 $\pm$ 0.000079	55.52 $\pm$ 9.48	4.67 $\pm$ 0.74							C
J123745.73+621456.55	189.440523 $\pm$ 0.000009	62.249043 $\pm$ 0.000009	20.30 $\pm$ 2.14	15.49 $\pm$ 1.03	16.86 $\pm$ 1.47	13.08 $\pm$ 0.71	16.86 $\pm$ 1.47	154 $\pm$ 27	75 $\pm$ 37	1	C
J123745.33+622023.75	189.438871 $\pm$ 0.000023	62.339931 $\pm$ 0.000053	19.11 $\pm$ 5.38	8.89 $\pm$ 1.78	10.04 $\pm$ 1.90	10.95 $\pm$ 1.15	10.48 $\pm$ 0.65	0 $\pm$ 88	0 $\pm$ 18	0	C
J123744.68+621218.76	189.436167 $\pm$ 0.000047	62.205212 $\pm$ 0.000044	57.75 $\pm$ 12.35	12.72 $\pm$ 2.27	20.58 $\pm$ 5.99	7.53 $\pm$ 1.66	20.58 $\pm$ 5.99	383 $\pm$ 125	203 $\pm$ 88	1	C
$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$
J123642.22+621545.48	189.175896 $\pm$ 0.000003	62.262634 $\pm$ 0.000003	54.77 $\pm$ 2.14	45.02 $\pm$ 1.08	47.62 $\pm$ 1.89						M
J123629.01+621045.59	189.120878 $\pm$ 0.000026	62.179330 $\pm$ 0.000054	26.55 $\pm$ 4.71	7.29 $\pm$ 1.04							M
J123612.46+621140.48	189.051920 $\pm$ 0.000024	62.194578 $\pm$ 0.000027	26.28 $\pm$ 4.17	8.72 $\pm$ 1.07	9.63 $\pm$ 1.99						M
J123531.57+621117.51	188.881551 $\pm$ 0.000023	62.188198 $\pm$ 0.000034	19.03 $\pm$ 4.53	10.70 $\pm$ 1.72	12.24 $\pm$ 2.71						M

References. — ^aThe rms position uncertainties are given by $\phi/[(2\ln 2)^{1/2}\times\rm SNR$] (Condon et al., 1998). ^b$S_{\rm I}$ – Integrated flux density. ^c$S_{\rm P}$ – Peak brightness. ^d$S_{\star}$ – Best estimate for the integrated flux density (see Section 3.3). ^eIf sources are confidently resolved along the major axis in the $0\hbox{$.\!\!^{\prime\prime}$}22$ resolution mosaic ($\phi_{\rm M}-\theta_{\rm 1/2}\geq 2\sigma_{\rm\phi_{\rm M}}$), $\rm R=1$. Else, if $\phi_{\rm M}-\theta_{\rm 1/2}<2\sigma_{\rm\phi_{\rm M}}$, $\rm R=0$. ^fSource type. $\rm T=C$ for compact, single radio sources in the $1\hbox{$.\!\!^{\prime\prime}$}0$ resolution mosaic and $\rm T=M$ for multi-component, extended and complex radio sources in the $1\hbox{$.\!\!^{\prime\prime}$}0$ or $0\hbox{$.\!\!^{\prime\prime}$}22$ resolution mosaic (see Tables A and A).

Table 4: Flux densities and positions of the radio sources in the low-resolution ($1\hbox{$.\!\!^{\prime\prime}$}0$ resolution) mosaic that are part of the 12 multi-component radio sources in our master catalog.

Name	RA1.0${}^{\tablenotemark{a}}$	DEC1.0${}^{\tablenotemark{a}}$	$S_{\rm I}^{\rm 1.0}$ ${}^{\tablenotemark{b}}$	$S_{\rm P}^{\rm 1.0}$ ${}^{\tablenotemark{c}}$
	(degs)	(degs)	$(\rm\mu Jy)$	$(\rm\mu Jy\,beam^{-1})$
J123820.47+621828.25	189.585301 $\pm$ 0.000006	62.307847 $\pm$ 0.000006	171.20 $\pm$ 7.71	86.55 $\pm$ 2.10
J123725.92+621128.38	189.359143 $\pm$ 0.000017	62.191054 $\pm$ 0.000020	92.86 $\pm$ 7.49	18.83 $\pm$ 1.28
	189.358146 $\pm$ 0.000002	62.191315 $\pm$ 0.000003	172.01 $\pm$ 3.83	99.40 $\pm$ 1.32
	189.356857 $\pm$ 0.000010	62.191148 $\pm$ 0.000009	106.40 $\pm$ 5.81	31.50 $\pm$ 1.36
J123717.90+621855.63	189.324579 $\pm$ 0.000020	62.315452 $\pm$ 0.000025	13.45 $\pm$ 2.77	7.83 $\pm$ 1.09
J123711.96+621325.92	189.299983 $\pm$ 0.000017	62.223809 $\pm$ 0.000014	9.50 $\pm$ 1.95	8.39 $\pm$ 1.03
	189.299380 $\pm$ 0.000046	62.224119 $\pm$ 0.000051	5.92 $\pm$ 2.56	3.72 $\pm$ 1.05
J123711.31+621330.93	189.297126 $\pm$ 0.000029	62.225258 $\pm$ 0.000012	26.74 $\pm$ 3.58	10.54 $\pm$ 1.04
J123707.99+621121.65	189.283196 $\pm$ 0.000024	62.189335 $\pm$ 0.000019	14.78 $\pm$ 2.80	8.00 $\pm$ 1.05
	189.284158 $\pm$ 0.000042	62.189409 $\pm$ 0.000059	5.26 $\pm$ 2.42	3.37 $\pm$ 1.01
J123645.81+620754.29	189.191129 $\pm$ 0.000023	62.131740 $\pm$ 0.000011	15.35 $\pm$ 2.54	9.77 $\pm$ 1.06
	189.190348 $\pm$ 0.000016	62.131746 $\pm$ 0.000026	6.85 $\pm$ 1.88	6.29 $\pm$ 1.01
J123644.40+621133.19	189.184950 $\pm$ 0.000001	62.192536 $\pm$ 0.000001	465.91 $\pm$ 3.21	402.52 $\pm$ 1.02
	189.184883 $\pm$ 0.000012	62.191452 $\pm$ 0.000045	11.98 $\pm$ 2.70	6.83 $\pm$ 1.01
J123642.22+621545.48	189.175896 $\pm$ 0.000003	62.262634 $\pm$ 0.000003	54.77 $\pm$ 2.14	45.02 $\pm$ 1.08
J123629.01+621045.59	189.120427 $\pm$ 0.000017	62.179271 $\pm$ 0.000026	10.71 $\pm$ 2.40	7.26 $\pm$ 1.05
	189.121352 $\pm$ 0.000037	62.179412 $\pm$ 0.000049	17.31 $\pm$ 4.14	5.45 $\pm$ 1.02
J123612.46+621140.48	189.051920 $\pm$ 0.000024	62.194578 $\pm$ 0.000027	26.28 $\pm$ 4.17	8.72 $\pm$ 1.07
J123531.57+621117.51	188.881551 $\pm$ 0.000023	62.188198 $\pm$ 0.000034	19.03 $\pm$ 4.53	10.70 $\pm$ 1.72

References. — ^aThe rms position uncertainties are given by $\phi/[(2\ln 2)^{1/2}\times\rm SNR$] (Condon et al., 1998). ^b$S_{\rm I}$ – Integrated flux density. ^c$S_{\rm P}$ – Peak brightness.

Table 5: Flux densities, sizes, and positions of the radio sources in the $0\hbox{$.\!\!^{\prime\prime}$}22$ resolution mosaic that are part of the 11 multi-component radio sources in our master catalog that have a counterpart in the high-resolution map.

Name	RA0.22${}^{\tablenotemark{a}}$	DEC0.22${}^{\tablenotemark{a}}$	$S_{\rm I}^{\rm 0.22}$ ${}^{\tablenotemark{b}}$	$S_{\rm P}^{\rm 0.22}$ ${}^{\tablenotemark{c}}$	$S_{*}^{\rm 0.22}$ ${}^{\tablenotemark{d}}$	$\theta_{\rm M}$	$\theta_{\rm m}$	R${}^{\tablenotemark{e}}$
	(degs)	(degs)	$(\rm\mu Jy)$	$(\rm\mu Jy\,beam^{-1})$	$(\rm\mu Jy)$	(mas)	(mas)
J123820.47+621828.25	189.585226 $\pm$ 0.000012	62.307875 $\pm$ 0.000011	196.28 $\pm$ 16.97	30.03 $\pm$ 1.45	196.28 $\pm$ 16.97	1169 $\pm$ 140	167 $\pm$ 33	1
	189.584946 $\pm$ 0.000004	62.308504 $\pm$ 0.000005	12.04 $\pm$ 3.03	9.80 $\pm$ 1.52	10.86 $\pm$ 0.69	132 $\pm$ 81	73 $\pm$ 109	0
J123725.92+621128.38	189.358115 $\pm$ 0.000001	62.191306 $\pm$ 0.000001	38.35 $\pm$ 1.79	32.74 $\pm$ 0.93	38.35 $\pm$ 1.79	106 $\pm$ 16	74 $\pm$ 20	1
	189.358420 $\pm$ 0.000001	62.191346 $\pm$ 0.000002	71.16 $\pm$ 3.34	31.70 $\pm$ 0.90	71.16 $\pm$ 3.34	311 $\pm$ 18	125 $\pm$ 15	1
	189.357841 $\pm$ 0.000001	62.191286 $\pm$ 0.000001	55.16 $\pm$ 2.84	26.40 $\pm$ 0.97	55.16 $\pm$ 2.84	251 $\pm$ 17	208 $\pm$ 16	1
J123717.90+621855.63	189.324547 $\pm$ 0.000006	62.315484 $\pm$ 0.000004	4.29 $\pm$ 1.30	3.98 $\pm$ 0.70	4.13 $\pm$ 0.67	0 $\pm$ 97	0 $\pm$ 29	0
	189.324573 $\pm$ 0.000005	62.315406 $\pm$ 0.000010	5.36 $\pm$ 1.78	3.38 $\pm$ 0.73	4.26 $\pm$ 0.70	272 $\pm$ 117	0 $\pm$ 40	0
J123711.96+621325.92	189.299945 $\pm$ 0.000003	62.223830 $\pm$ 0.000004	8.90 $\pm$ 1.55	6.78 $\pm$ 0.74	7.77 $\pm$ 0.72	183 $\pm$ 54	33 $\pm$ 148	0
J123711.31+621330.93	189.297436 $\pm$ 0.000004	62.225281 $\pm$ 0.000005	6.94 $\pm$ 1.61	5.13 $\pm$ 0.75	5.97 $\pm$ 0.73	163 $\pm$ 72	94 $\pm$ 84	0
	189.296870 $\pm$ 0.000002	62.225240 $\pm$ 0.000002	8.02 $\pm$ 1.16	8.63 $\pm$ 0.70	8.32 $\pm$ 0.67	0 $\pm$ 80	0 $\pm$ 15	0
J123707.99+621121.65	189.283196 $\pm$ 0.000008	62.189312 $\pm$ 0.000006	8.26 $\pm$ 2.11	3.94 $\pm$ 0.72	8.26 $\pm$ 2.11	293 $\pm$ 91	169 $\pm$ 76	1
J123645.81+620754.29	189.191093 $\pm$ 0.000015	62.131730 $\pm$ 0.000007	12.33 $\pm$ 3.17	3.40 $\pm$ 0.70	12.33 $\pm$ 3.17	502 $\pm$ 138	233 $\pm$ 80	1
	189.190412 $\pm$ 0.000009	62.131725 $\pm$ 0.000011	12.40 $\pm$ 2.88	3.93 $\pm$ 0.71	12.40 $\pm$ 2.88	466 $\pm$ 118	198 $\pm$ 70	1
J123644.40+621133.19	189.184941 $\pm$ 0.000001	62.192538 $\pm$ 0.000001	414.04 $\pm$ 2.82	339.97 $\pm$ 0.70	414.04 $\pm$ 2.82	250 $\pm$ 1	66 $\pm$ 2	1
J123642.22+621545.48	189.175882 $\pm$ 0.000001	62.262644 $\pm$ 0.000001	44.68 $\pm$ 1.33	41.63 $\pm$ 0.73	44.68 $\pm$ 1.33	78 $\pm$ 13	34 $\pm$ 25	1
	189.175916 $\pm$ 0.000010	62.262546 $\pm$ 0.000006	2.89 $\pm$ 1.34	2.66 $\pm$ 0.72	2.77 $\pm$ 0.69	0 $\pm$ 148	0 $\pm$ 47	0
J123612.46+621140.48	189.051956 $\pm$ 0.000004	62.194610 $\pm$ 0.000004	6.69 $\pm$ 1.49	5.37 $\pm$ 0.74	5.99 $\pm$ 0.71	132 $\pm$ 72	82 $\pm$ 89	0
	189.051860 $\pm$ 0.000006	62.194679 $\pm$ 0.000013	2.95 $\pm$ 1.33	2.63 $\pm$ 0.67	2.79 $\pm$ 0.64	0 $\pm$ 155	0 $\pm$ 31	0
J123531.57+621117.51	188.881592 $\pm$ 0.000004	62.188180 $\pm$ 0.000005	5.60 $\pm$ 1.82	6.13 $\pm$ 1.11	5.86 $\pm$ 0.65	0 $\pm$ 312	0 $\pm$ 34	0
	188.881497 $\pm$ 0.000006	62.188251 $\pm$ 0.000005	6.64 $\pm$ 2.02	6.46 $\pm$ 1.10	6.55 $\pm$ 0.65	175 $\pm$ 97	0 $\pm$ 24	0

References. — ^aThe rms position uncertainties are given by $\phi/[(2\ln 2)^{1/2}\times\rm SNR$] (Condon et al., 1998). ^b$S_{\rm I}$ – Integrated flux density. ^c$S_{\rm P}$ – Peak brightness. ^d$S_{\star}$ – Best estimate for the integrated flux density (see Section 3.3). ^eIf sources are confidently resolved along the major axis in the $0\hbox{$.\!\!^{\prime\prime}$}22$ resolution mosaic ($\phi_{\rm M}-\theta_{\rm 1/2}\geq 2\sigma_{\rm\phi_{\rm M}}$), $\rm R=1$. Else, if $\phi_{\rm M}-\theta_{\rm 1/2}<2\sigma_{\rm\phi_{\rm M}}$, $\rm R=0$.