Halfway to the Peak: The JWST MIRI 5.6 micron number counts and source population

Abstract

We present an analysis of 8 JWST Mid-Infrared Instrument 5.6 micron images with $5\,\sigma$ depths of $\approx$0.1$\mu$Jy. We detect 2854 sources within our combined area of 18.4 sq.arcmin – a $>4\times$ increase in source density over earlier IRAC channel 3 data. We compute the MIRI 5.6 $\mu$m number counts including an analysis of the field-to-field variation. Relative to earlier published MIRI 5.6 $\mu$m counts, our counts have a more pronounced knee, at roughly 2 $\mu$Jy. The location and amplitude of the counts at the knee are consistent with the Cowley et al. (2018) model predictions, although these models tend to overpredict the counts below the knee. In areas of overlap, 84 % of the MIRI sources have a counterpart in the COSMOS2020 catalog. These MIRI sources have redshifts that are mostly in the $z\sim 0.5-2$, with a tail out to $z\sim 5$. They are predominantly moderate to low stellar masses ($10^{8}-10^{10}$M_⊙) main sequence star-forming galaxies suggesting that with $\approx$2 hr exposures, MIRI can reach well below $M^{*}$ at cosmic noon and reach higher mass systems out to $z\sim 5$. Nearly 70% of the COSMOS2020 sources in areas of overlap now have a data point at 5.6 $\mu$m (rest-frame near-IR at cosmic noon) which allows for more accurate stellar population parameter estimates. Finally, we discover 31 MIRI-bright sources not in COSMOS2020. A cross-match with IRAC channel 1 suggests that 10-20 % of these are likely lower mass (M${}_{*}\approx 10^{9}$M_⊙), $z\sim 1$ dusty galaxies. The rest (80–90%) are consistent with more massive, but still very dusty galaxies at $z>3$.

revtex4-1.clsrevtex4-2.cls

1 Introduction

Extragalactic astronomy aims to study how galaxies form and evolve across time. To do so, we build multiwavelength galaxy surveys which allow us to compare the data with spectral energy distribution (SED) models to infer redshifts and stellar population parameters (e.g. Brammer et al., 2008; Weaver et al., 2022a; Wang et al., 2024). Mid-IR data are particularly critical to this effort since they probe the rest-frame near-IR in the cosmic noon epoch and earlier which allows for much more accurate stellar masses and star formation rates (SFRs) (e.g. Elsner et al., 2008; Stefanon et al., 2015; Muzzin2009l; Martis et al., 2023; La Torre et al., 2024). Indeed, over the past two decades Spitzer IRAC data – especially its first two channels, available even in the warm mission (Lacy et al., 2021; Annunziatella et al., 2023) have been a critical component of these multiwavelength surveys. With IRAC channel 1 and 2, we reached stellar masses of $\approx 10^{9.5}$M$\odot$ at cosmic noon, i.e. well below the knee of the SMF (Elsner et al., 2008; Madau & Dickinson, 2014), but only detected the most massive galaxies ($log(M_{*})>11$) at $z>4$ (e.g. Stefanon et al., 2015). The IRAC channel 3, at 5.8 $\mu$m, which critically covers the rest-frame of the stellar 1.6 $\mu$m bump at $z\sim 2-3$ was relatively underutilized because it was not available during the extended warm mission and was the least sensitive IRAC channel even during the cold mission.

The JWST (Gardner et al., 2006), which launched on December 25, 2021, is already revolutionizing infrared astronomy with its unprecedented performance (Rigby et al., 2023). In particular, the Mid-InfraRed Instrument (MIRI; Rieke et al., 2015) allows for significantly greater sensitivity and angular resolution relative to the Spitzer Space Telescope. Recently published JWST/MIRI number counts at 7.7 $\mu$m, 10 $\mu$m and 15 $\mu$m bands (Ling et al., 2022; Wu et al., 2023; Kirkpatrick et al., 2023; Stone et al., 2024) show dramatic improvement in the depth reached, even with moderate exposure times, relative to prior measurements from Spitzer and ISO. These counts have already been used in constraining galaxy and black hole evolution models (Kim et al., 2024). Number counts at the shortest MIRI wavelength (5.6 $\mu$m) are much more scarse (Yang et al., 2023; Stone et al., 2024), but they are critical as this band samples the stellar 1.6 $\mu$m bump at $z\sim 2-3$ and thus is critical in testing our models of the galaxy population at cosmic noon.

In addition, sampling the rest-frame near-IR at cosmic noon, the 5.6 $\mu$m band is much less sensitive to the effects of dust obscuration even than traditionally ‘dust clean’ bands such as the $K$ band, which is rest-frame $r$-band at the same redshifts. This insensitivity to dust obcuration is important as in the more than two decades since the discovery of the Cosmic Infrared Background (e.g. Puget et al., 1996) it has become abundantly clear that the bulk of star-formation activity at cosmic noon and beyond takes place in dust obscured galaxies (see Casey et al. 2014 for a review) which in their extreme are missed in UV/optical surveys (see e.g. Hughes et al., 1998; Sajina et al., 2006; Pope et al., 2008; Zavala et al., 2021, among many others). Recently, JWST data have further highlighted this by finding that even deep HST data miss the reddest/most dust obscured parts of the galaxy population (e.g. Labbe et al., 2023; Barrufet et al., 2023; Williams et al., 2023).

In this paper, we use deep MIRI 5.6 $\mu$m images to provide a first look at the properties of the 5.6 $\mu$m number counts and source population based on images obtained with nearly 2 hr exposure times. This wavelength corresponds roughly to the rest-frame $H/K$-band at cosmic noon ($z\sim 1-3$), thus probing primarily stellar mass at these critical redshifts. This paper is structured as follows. In Section 2 we present the MIRI 5.6 $\mu$m imaging data and data reduction. In Section 3 we present the source detection and photometry and the verification of the latter. In Section 4, we show the key results of our work. These include the 5.6 $\mu$m number counts (Section 4.2); the redshift distribution and other properties of the MIRI 5.6 $\mu$m sources with counterparts in COSMOS2020 (Section 4.3); and lastly MIRI 5.6 $\mu$m sources without counterparts in COSMOS2020 suggesting very red colors (Section 4.4). Throughout we adopt the AB magnitude system (Oke & Gunn, 1974). We adopt the cosmology model from the COSMOS2020 catalog (Weaver et al., 2022a), which is a standard $\Lambda$CDM cosmology with H₀ = 70 km/s/Mpc, $\Omega_{M,0}$ = 0.3, and $\Omega_{\Lambda,0}$ = 0.7.

2 Data

2.1 MIRI imaging data

Our data come from JWST/MIRI imaging with the 5.6 $\mu$m filter obtained as part of the GO1 program “Halfway to the Peak: A Bridge Program To Map Coeval Star Formation and Supermassive Black Hole Growth” (PIs Pope, Sajina, Yan; PID 1762; Pope et al., 2021). This program observed 8 targets with the MIRI/MRS instrument: details on the targets and MRS spectra analysis are provided in Young et al. (2023). Simultaneous with the MIRI/MRS observations, we obtained MIRI 5.6 $\mu$m imaging in fields adjacent to the MRS targets. This parallel observing means our imaging fields are effectively blank fields and therefore ideal for statistical studies such as number counts. In total, we have six fields in the First Look Survey (FLS) field and two in the COSMOS field. The FLS fields were observed in July/August 2022 and the COSMOS fields were observed in December 2022.

Field name	RA	Dec	Exp. time
FLS1	17:12:28.50	+58:59:30.12	1.85hrs
FLS2	17:24:46.96	+59:15:24.01	1.85hrs
FLS3	17:21:07.16	+58:45:39.67	1.85hrs
FLS4	17:22:49.98	+59:40:32.35	1.85hrs
FLS5	17:19:12.87	+59:28:53.55	1.85hrs
FLS6	17:13:13.14	+58:55:22.51	1.85hrs
COSMOS1	10:01:14.98	+2:24:36.36	1.85hrs
COSMOS2	10:01:42.53	+2:47:26.64	1.85hrs

Table 1: The centers of the MIRI imaging fields used in this paper.

We obtain data over the MIRI imager with an FOV of $74\times$113^′′ as well as the smaller MIRI coronograph with a FOV of $24\times$24^′′ (as seen in Figure 1). For our eight pointings, this adds up to a total area of just under 19.4 sq.arcmin. The effective area we use is a bit more restricted after masking out noisy edges and bright stars. This is discussed in Section 3.2.

Figure 1: An example of our MIRI images (for the FLS6 field) highlighting the number of sources clearly detected within the imager’s FOV of only 2.3 sq.arcmin. The angular resolution is also illustrated in the zoomed-in insets, both measuring 9$\times$9^′′, where we see both resolved morphology as illustrated by the detected spiral structure (top inset) as well as resolved closely spaced sources (that would have been blended in the Spitzer IRAC images, middle inset). The 2^′′ scale bars are comparable to the IRAC channel 3 FWHM which is 1.88^′′.

Figure 2: Example images produced with the standard pipeline, and following two different background removal methods, as described in Section 2.2. The image median (in counts/second/pixel) are given in the bottom-left of each panel. Left: the FLS6 image from the standard pipeline. Middle: the FLS6 image where we have subtracted a median background constructed from all six FLS fields, with the sources masked out. Right: the FLS6 image with super background subtraction (Pérez-González et al., 2024; Alberts et al., 2024).

2.2 Data reduction

We first obtained our MIRI images from the archive already reduced using the JWST operations pipeline build 10.1 (science calibration pipeline version 1.12.5; CRDS version 11.17.6, context jwst_1193.pmap). Figure 2 left shows an example of one of our images. Within these pipeline-reduced images we noted both the presence of an overall background level (the median pixel value for the FLS images was 1.23, for the COSMOS images it was 1.52). This difference in overall background level is unsurprising given that FLS is a northern field, whereas COSMOS is equatorial and as such sees a higher level of Zodiacal light, see example discussion in Sanders et al. (2007). We also noted some clear striping especially toward the bottom quarter of each image.

We first attempted a simple background subtraction motivated by the fact that the stripping patterns were fairly consistent among the images. We thus constructed a median background image where we used the segmentation map provided with pipeline-processed images to mask out bright sources in the fields. To avoid any remaining source signal, we added a 5 pixel “buffer” around each source within the segmentation maps. Due to the difference in overall pixel level between the images, we rescaled each source-subtracted image to a median of one. Figure 2 middle shows an example of our image after this basic median background subtraction when we only use the FLS images for constructing the background. It is clear that this method works well in removing the vast majority of the striping though some residual structure remains due to second order differences between images. In addition, since we only have two COSMOS images we could not construct a COSMOS only median using this method. We attempted a combined median background in FLS and COSMOS; however, this led to the COSMOS fields (which have less weight in this median with 2/8 images only) to have higher noise levels. Specifically, using a single median we found 5 $\sigma$ levels of $\approx$0.095 $\mu$Jy in FLS and $\approx$0.135 $\mu$Jy in COSMOS. When using only the FLS images for a median (as shown in Figure 2 middle), we find 5 $\sigma$ levels of $\approx$0.085 $\mu$Jy. While we accounted for the differences in overall background level, these results clearly indicate that we have some time variable aspects to the background which make it sub-optimal to construct such median background from observations about 6 months apart.

To address these issues, we then processed our images through the super-background subtraction procedure developed by the SMILES (PID 1207; PI G. Rieke; Lyu et al. (2024)) team (Pérez-González et al., 2024; Alberts et al., 2024). This background subtraction starts with the stage 2 images and homogenizes the background across all images in the program taking into account the time varying behavior of the background which is particularly prominent at shorter wavelengths as here. The full details are described in Alberts et al. (2024). Figure 2 right shows the result of this super background subraction. It is clear that the residual structure in our previous simple background-subtracted image is now largely gone. In addition, this method allows for COSMOS only median since it no longer uses the fully processed images (of which there are only two), but the stage 2 images of which we have 2$\times$6 dithers which is sufficient for a median background. As we show in Section 2.4 below, we also no longer find a difference in the obtained depth in FLS and COSMOS. The depth we find in FLS relative to the simple background subtraction using the FLS images alone however are comparable. This suggests that depending on the specific science applications and data available, a simple procedure as described above may be sufficient, but it should only be applied to images obtained close in time to avoid time-varying background effects.

2.3 Noise properties

To explore the noise properties of our images post background subtraction, we constructed a noise map where each pixel value equals the standard deviation in the image pixel values in a 3$\times$3 pixel grid around that pixel. We performed this calculation on each image with two masks applied. First, we masked the pixels belonging to sources by applying the segmentation map. We then created a mask that removes the image edges and the residual strips, as we consider them to be sources of additional noise. This second mask is the same for each image. We perform our noise calculation on the science images after the segmentation map and noise/stripe mask have been applied. The resulting noise map for the FLS1 field is shown in Figure 3. We compared our produced noise map to the one provided with the standard pipeline products (which is derived from the per pixel coverage). The typical pixel values of these calculated noise maps are 60-70% of the MAST-provided error maps.

Figure 3 shows a wavy pattern across the whole image. We found this pattern in all our images, even when varying the number of pixels used in calculating the standard deviation, and also when running the same procedure before or after background subtraction. We also found this pattern with both the simple and super background subtracted images. We did not find this pattern when running the same procedure on one of the publicly-available CEERS images in the same filter (Finkelstein et al., 2022). The CEERS images have lower exposure times (0.82 vs. 1.85 hrs), and slightly smaller pixel scales (0.09^′′ vs 0.11^′′). At present, we consider this an unexplained instrumental effect that is consistent with the description of “tree rings” that have previously been observed in MIRI imaging. This pattern is at a very low level, with the amplitude of the peaks and troughs translating to $\approx$ 0.5 nJy/pixel, as can also be seen in the scale color bar.

2.4 Image depth

We estimate the image depth by placing random apertures in source-free parts of the background-subtracted images. The standard deviation of the Gaussian distribution of the resultant empty aperture fluxes is the $1\,\sigma$ depth of the image. The above method is a simplification of the method used in Annunziatella et al. (2023), since we have fairly uniform coverage unlike that in Annunziatella et al. (2023) so noise-weighting is not critical. In Section 3.1 we describe our SExtractor photometry. For our random apertures we use a diameter equal to the median of the SExtractor Kron radii, where a Kron radius is a “reduced pseudo-radius” that defines a circular aperture containing $>90\%$ of the flux of the source. We choose Kron apertures as they are meant to capture most of the light of our sources, which are typically spatially extended. The median Kron radius in our catalog is $\sim$0.385^′′, so we adopt a diameter of 0.77^′′ for our apertures. For comparison, the FWHM of the MIRI PSF at 5.6 $\mu$m is 0.207^′′.

We placed 500 empty apertures per field for a total of 4000 apertures. For each image, we then fitted a Gaussian to the resulting histogram of aperture flux. We found that the standard deviations of these fitted Gaussians were consistent between all the images, translating to a $5\,\sigma$ depth of 0.1 $\mu$Jy.

As a cross-check, the median flux uncertainties for our sources based on the SExtractor photometry (Section 3.1) translate to a 5 $\sigma$ uncertainty of 0.081 $\mu$Jy. This is reasonably consistent with our image depth analysis above. Note that the random apertures method is sensitive to the adopted aperture sizes and slightly smaller apertures would bring the two values into even closer agreement. As an additional cross-check, the image depths quoted by CEERS (Yang et al., 2023), rescaled to our exposure time (since depth scales as $\sqrt{t}$) translate to $\approx$0.08 $\mu$Jy for our images. To be conservative, we adopt $\approx$0.1 $\mu$Jy as our nominal 5 $\sigma$ image depth.

Figure 3: Noise map for the FLS1 field. Produced by computing the standard deviation of the FLS1 image inside a $3\times 3$ pixel square centered on each pixel. The science image used in producing the noise map was first multiplied by the field’s segmentation map, removing the sources. We then applied a mask, shown here in orange, that removes edges and residual stripes as additional sources of noise. After both steps, we produced the noise map.

2.5 Photometric zero-point

Pre-launch estimates of the MIRI imager’s performance already predicted that its point spread function (PSF) has a ‘cruciform’ artifact that draws luminosity from the core into the wings of the Gaussian. In-flight measurements have shown that this artifact is prominent at shorter wavelengths like ours (Rigby et al., 2023), but is negligible at $>$10 $\mu$m (Gáspár et al., 2021). The extent of this effect was only noted once in-flight assessments had been made. Since our data were obtained very early after the start of science operations, the FITS header keyword PHOTUJA2 provided for converting from the native units to physical photometry units was incorrect. The fiducial keyword translated to a photometric zero-point of 26.121 mag which resulted in us seeing a significant systematic offset in the photometry measured by MIRI F560W as compared to the Spitzer IRAC Ch3 data available for our fields.

This issue was also noted by the CEERS team and led them to calculate a photometric zero-point for the 5.6 $\mu$m images of 25.701 in AB magnitudes (Papovich et al., 2023). However, the CEERS images have a pixel scale of 0.09$\arcsec$ (Papovich et al., 2023) whereas ours is 0.11$\arcsec$. We re-scale the CEERS calculated zero-point as in $zp=25.701-5\log_{10}(0.11/0.09)=25.265$. This is the value we adopt in this paper which leads to much better agreement between the MIRI 5.6 $\mu$m and IRAC channel 3 photometry for cross-matched sources (see Section 3.3).

3 Data analysis

3.1 SExtractor photometry

We used Source Extractor (Bertin, E. & Arnouts, S., 1996) version 2.28.0 for source detection and photometry. Table 2 lists our SExtractor parameter settings. The SEEING$\_$FWHM keyword adopted here corresponds to the FWHM of the in-flight measured MIRI 5.6 $\mu$m PSF (see section 4.1 for details). As our fiducial flux measurements, we adopt values produced with Kron aperture photometry, labeled by the FLUX_AUTO keyword by SExtractor. The photometric uncertainties are natively provided by SE and take into account the variance of the pixel values inside the aperture, the background value, and the effective detector gain.

Extraction
DETECT_MINAREA	9.0
DETECT_THRESH	1.5
ANALYSIS_THRESH	0.4
FILTER_NAME	gauss_3.0_5x5.conv
Photometry
PHOT_AUTOPARAMS	2.5, 3.5
PHOT_FLUXFRAC	0.5
MAG_ZEROPOINT	25.265
GAIN	36630.528
PIXEL_SCALE	0.11091469
Background
BACK_TYPE	AUTO
BACK_SIZE	16
BACK_FILTERSIZE	3
BACKPHOTO_TYPE	LOCAL
Deblending
DEBLEND_NTHRESH	32
DEBLEND_MINCONT	0.003
Star/Galaxy Separation
SEEING_FWHM	0.207

Table 2: SExtractor parameter settings for the photometry.

3.2 Masking

We further clean the raw SExtractor catalogs by masking out bright stars as well as any particularly noisy parts of the images. To do so we created masks for each image. These masks include two components. One excludes obviously diffracted sources. The other avoids all edges (assuming an edge width of 5 pixels), where edges include the blocked regions of the Lyot coronograph. We also manually remove the noisy residual stripes near the bottom of the images (see right-hand panel in Figure 2). Figure 3 shows as an example the mask for the FLS1 field.

We also apply a quality cut on our catalogs, only keeping sources with a signal-to-noise ratio greater than 5. The combination of the masking and the quality cut reduces the raw SExtractor output of 4,089 sources to a total of 2,854 “reliable” sources (“reliable” meaning falling outside our defined masks and having SNR$\geq 5$) among all our eight fields. Per field, we detect the most in field FLS3 (409 sources) and the least in field FLS5 (310 sources).

Figure 4: Comparison of the position centroids of the 582 matched sources between MIRI and the COSMOS2020. The histograms (top and right) represent the offsets in both fields combined, with the mean absolute deviation (MAD, calculated as the mean of $|\Delta\mathrm{RA}|$ and $|\Delta\mathrm{Dec}|$, respectively) labeled.

3.3 Astrometric & photometric verification

In order to perform astrometric and photometric verifications, we used the SExtractor detected sources in our two COSMOS fields and cross-matched them with the COSMOS2020 CLASSIC catalog (Weaver et al., 2022a). The COSMOS2020 source detection is based on a weighted combined $izYJHK$ image (Weaver et al., 2022a). We used a 1^′′ matching radius and found a total of 582 COSMOS2020 cross-matches, 175 of which also have IRAC channel 3 detections. Overall, our MIRI data increase the fraction of COSMOS2020 with 5.6/5.8 $\mu$m data by roughly 3 times. Note that our MIRI fields are within the deep part of the COSMOS2020 catalog and not near the edges and the COSMOS2020 positions rely on Gaia astrometric solutions (Weaver et al., 2022a). Figure 4 shows the comparison between the MIRI and COSMOS2020 source positions of the cross-matched sources. The offsets we find are small and well within the MIRI PSF at this wavelength. Therefore we do not explicitly apply astrometric corrections.

Figure 5 shows the comparison between the MIRI 5.6 $\mu$m and the available IRAC channel 3 fluxes. As described in section 3.1, we use Kron aperture photometry for the 5.6 $\mu$m fluxes. For the IRAC fluxes, we used the aperture corrected 2^′′ aperture photometry from Weaver et al. (2022b). We overlay the 1:1 line. Overall, we find good agreement, although we have a median offset of 0.22 magnitudes.

Figure 5: Photometry cross-check for the matched sources with available IRAC channel 3 measurements applying a $\mathrm{SNR}>1$ cut on the IRAC ch3 data, leaving us with 119 sources. The dashed diagonal line is the 1:1 line. The vertical dashed line represents the $3\sigma$ IRAC channel 3 limit (Weaver et al., 2022a). The quoted median offset is computed for the SNR$>3$ sources.

3.4 Estimating fraction of fake sources

The cross-match with COSMOS2020 also gives us an upper limit on the fraction of the SExtractor detected sources that might be fake. Such sources are expected due to some of the residual artifacts in the images, that may have been missed in the masks described in Section 3.2. Overall we found COSMOS2020 counterparts for 582 of the 690 MIRI sources in our two COSMOS fields (84 %). The 108 unmatched sources ($16\,\%$), represent an upper limit on the fake source fraction since it is expected that we have some red MIRI sources that are not detected in COSMOS2020. Figure 6 shows the 5.6 $\mu$m flux distributions of the COSMOS2020 undetected sources vs. the matched sources. Unsurprisingly, the unmatched sources are generally fainter than the overall MIRI source population. Indeed, visual examination of the 108 unmatched sources showed only 31 are unambiguous MIRI detections. Removing these sources which are clean MIRI detections despite being unmatched in COSMOS2020 leaves us with 77 likely fake sources ($\sim 10\%$ of the total). Many of these apparently fake sources are driven by single bright pixels that were missed in the standard pipeline reduction, others appear to be residuals of the striping pattern seen in Figure 2. In Figure 6 we show separately the flux distributions of the unmatched sources judged to be real vs. those that are potentially fake. Note that the higher Section 4.4 explores in more detail the nature of the 31 reliable MIRI sources that are unmatched in COSMOS2020.

Figure 6: A comparison between the flux distributions of our MIRI sources that are matched to COSMOS2020 sources vs. those that are unmatched.

4 Results

4.1 Source completeness

In order to estimate the source detection completeness at different flux density levels, we take the usual Monte Carlo approach by injecting fake sources at different flux levels and noting the fraction of them that are recovered with the same SExtractor procedure as applied to real sources (e.g. Takagi et al., 2012). We use the code described in Shipley et al. (2018) for this step.

For the injected fake sources, we need a model of their on-the-sky distribution. Simply adopting the MIRI PSF, is not appropriate since the vast majority of our sources are spatially resolved. To create a more realistic source model, we stack 500 isolated and visually compact sources to construct a source model stamp. Figure 7 shows this stamp which measures 16 pixels $\times$ 16 pixels ($\approx$1.8${}^{\prime\prime}\times$1.8^′′). This source model stamp has a mean (along x- and y-axis) full width at half-maximum (FWHM) of 0.357^′′ which, as expected, is larger than the MIRI 5.6 $\mu$m PSF of 0.207^′′, as measured during JWST commissioning.

Figure 7: The $16\times 16$ pixel source model stamp used in the Monte Carlo simulations to determine source completeness. The FWHM of the best-fit Gaussian is 0.341^′′ in the x-direction and 0.330^′′ in the y-direction.

Completeness estimation starts with a segmentation map produced by SExtractor where all pixels with zero value correspond to no source present. We additionally set all pixel values to zero if they were masked, as described in Section 3.2, as well as setting to zero all sources consisting of 25 or fewer contiguous pixels. This second step was done in order to provide more sites suitable for injection, as our procedure avoids injecting in sites flagged by the segmentation map as having sources, though we deem this unnecessary for smaller sources. We then proceeded to inject fake sources across a range of flux densities ranging between 0.01 $\mu$Jy and 0.25 $\mu$Jy in 25 total increments. Sources are injected only on valid pixels not removed in the above steps. At each flux level, we injected 500 sources in the image. The sources are injected by using the source model scaled to the desired flux density. We then processed the newly-created images with the same SExtractor procedure used for the science images themselves. The ratio of recovered to injected sources at each flux density level constitutes the completeness at that level. This procedure was done 10 times per flux level per image, totaling 80 runs. Figure 8 shows the measured completeness curves for each run as well as the adopted completeness curve, which is the median between the 80 runs. As expected we find that the completeness is essentially 100% above our estimated 5 $\sigma$ image depth and drops rapidly below that. The 50% completeness is at $\approx$0.07 $\mu$Jy.

Figure 8: Source completeness for our images. The $x$-axis is trimmed to show detail. A sample of the computed completeness curves is faintly overlaid to show spread. The thin dotted line shows the $5\sigma$ depths of the images as determined by the empty apertures method, detailed in section 2.4.

Figure 9: Left The integral number counts of our individual 8 fields. For comparison we also show the combined counts with Poisson errors. We find that the typical field-to-field variation is a factor of roughly two – well about the simple Poisson errors on the combined counts. We also overlay the $5\,\sigma$ limit as a vertical dashed line.
Right Our integral number counts for our eight fields (orange symbols), with the Poisson and field-to-field variation uncertainties added in quadrature. We overlay for comparison both the published CEERS counts (Yang et al., 2023, blue crosses, include only Poisson error) and the SMILES counts (Stone et al., 2024, green circles, include only Poisson error), as well as model predictions MIRI F560W counts from Cowley et al. (2018, black curves) where the dot-dash curve is their fiducial model and the upper (dashed) one is the same but with evolving supernovae feedback (“EvolFB”). We also overlay the $5\,\sigma$ limit as a vertical dashed line.

4.2 Number counts

We calculate the MIRI 5.6 $\mu$m integral number counts by first counting the number of objects in each flux bin $i$ as $dN/dS_{i}$. The source completeness for each flux bin is $C_{i}$ and is as shown in Figure 8. The integral counts are given by:

$$N(>S)=\sum_{S_{i}=S}^{\infty}\frac{1}{A_{\mathrm{eff},i}}\frac{1}{C_{i}}\frac{dN}{dS_{i}}$$ (1)

We calculate the effective area, $A_{\mathrm{eff}}$, as the inverse of the Boolean mask used for source detection (Section 3.2) times the single pixel area. For all 8 fields combined, this adds up to 18.4 sq.arcmin. Note that our coverage is fairly uniform within these masks therefore we assume the effective area does not change with flux level. The calculated integral counts are given in Table 3 along with their Poisson uncertainties. To calculate the raw counts, i.e. the ones not corrected for source comleteness, we simply adopt $C_{i}\equiv 1$ in all flux bins in Equation 1. When computing our number counts, we use less stringent quality cut on our catalog, counting all sources with an SNR $\geq 1$. In Figures 9 we overlay the $5\sigma$ limit of $\approx$0.1 $\mu$Jy on our counts figures as reference.

Figure 9 left shows both our combined MIRI 5.6 $\mu$m integral number counts, with Poisson errors, and the counts estimated for each of our 8 fields individually. It is obvious that the field-to-field variation is significantly greater than the Poisson uncertainties on the combined counts. The shape of the counts is fairly consistent with nearly all fields showing a knee at around 6-8 $\mu$Jy.

Figure 9 right shows the combined integral MIRI 5.6 $\mu$m number counts where the uncertainties on the points represent the total uncertainties. The later are estimated by adding in quadrature the Poisson uncertainties with a field-to-field variation uncertainty which is taken to be half of the total spread between the individual field counts at each flux level. The combined counts and their Poisson and total uncertainties are given in Table 3.

For comparison, in Figure 9 we overlay the CEERS 5.6 $\mu$m counts from (Yang et al., 2023). The data for the CEERS number counts come from four nearby fields (separated by $<$15^′) with a total area of $\sim 9.5\,$sq.arcmin. The median $5\sigma$ depth for the CEERS images at 5.6 $\mu$m is 0.138 $\mu$Jy (Yang et al., 2023). We note that our counts tend to be higher than the CEERS ones, but this disagreement is only at $\approx$1 $\sigma$ level for most flux bins considering our total uncertainties that account for field-to-field variation. Note that the uncertainties plotted on the CEERS counts are purely Poisson and do not include field-to-field uncertainties. The gap between our counts and their only becomes significant at the lowest flux bins, below $\approx$0.3 $\mu$Jy; the widest part of the gap is $\sim 2.3\sigma$. We also overlay the just released MIRI 5.6 $\mu$m counts from the SMILES team (Stone et al., 2024). They are in good agreement with our counts at the bright end, but closer to CEERS at the faint end, though again within our 1 $\sigma$ uncertainties. We note that both CEERS and SMILES use combined F560W+F770W images for source detection, unlike this analysis which is done entirely with MIRI 5.6 $\mu$m images. To test the potential effect of differences in photometry, we ran our SExtractor setup on the public background-subtracted CEERS images and compared with their photometric catalog. We find excellent agreement, with the largest discrepancies per flux bin (e.g. for $<$ 1 $\mu$Jy fluxes) being that we find roughly 10% more sources. This is consistent with our estimated fake source fraction in Section 3.4. Such fake sources are less likely in the CEERS photometric catalog due to their including the F770W as well for source detection. This accounts for some of the offset we see at fainter fluxes, but $\approx$10% likely fake source fraction is far from sufficient to explain the gap which is $\approx$2 $\times$ at its maximum. Our Figure 9 left suggests this is primarily due to field-to-field variation instead.

In Figure 9 right we also overlay the model counts from Cowley et al. (2018). Their predictions are based on the GALFORM semi-analytic model for galaxy formation within a $\Lambda$CDM framework (Lacey et al., 2016). The predictions based on the core galaxy formation model from Lacey et al. (2016) are referred to as the ”Baseline model”. The model that includes supernovae feedback whose strength evolves with redshift is called the “EvolFB model”. The GALFORM model assumes the Maraston (2005) stellar population models and the dust radiative transfer comes from the GRASIL prescription of Silva et al. (1998). Our counts, in the brighter regime, are reasonably consistent with these model predictions with a slight preference for their evolving feedback model, although the difference between their models is well within our uncertainaties (including the field-to-field variation). In agreement with the CEERS and SMILES results, we find that the Cowley et al. (2018) models overestimate the counts for sources below $\approx$3 $\mu$Jy.

Table 3: The MIRI 5.6 $\mu$m number counts

log($S_{560W}$)	Raw $N(>S)$	Corrected $N(>S)$	Poisson error	Total error lower-upper
[$\mu$Jy]	[10-5deg${}^{-}2$]	[10-5deg${}^{-}2$]	[10-5deg${}^{-}2$]	[10-5deg${}^{-}2$]
0.05	5.49	6.82	0.10	1.37-1.40
0.07	5.40	6.01	0.10	0.96-1.05
0.09	5.01	5.10	0.10	0.63-0.60
0.12	4.41	4.41	0.09	0.54-0.55
0.16	3.79	3.79	0.09	0.58-0.55
0.22	3.23	3.23	0.08	0.49-0.53
0.31	2.65	2.65	0.07	0.50-0.46
0.42	2.25	2.25	0.07	0.58-0.28
0.58	1.87	1.87	0.06	0.58-0.30
0.80	1.54	1.54	0.05	0.45-0.30
1.10	1.32	1.32	0.05	0.38-0.21
1.51	1.12	1.12	0.05	0.35-0.23
2.07	0.94	0.94	0.04	0.29-0.24
2.85	0.80	0.80	0.04	0.28-0.22
3.92	0.68	0.68	0.04	0.26-0.14
5.39	0.54	0.54	0.03	0.22-0.17
7.41	0.43	0.43	0.03	0.22-0.13
10.18	0.31	0.31	0.02	0.13-0.06
14.00	0.21	0.21	0.02	0.08-0.08
19.25	0.13	0.13	0.02	0.07-0.08
26.47	0.07	0.07	0.01	0.04-0.07
36.39	0.05	0.05	0.009	0.03-0.05
50.04	0.02	0.02	0.006	0.0086-0.026
68.80	0.0058	0.0058	0.003	0.010-0.025
94.60	0.0019	0.0019	0.002	0.014-0.014

Figure 10: The COSMOS1 field with all COSMOS2020 Classic catalog (Weaver et al., 2022a) sources overlaid as the green circles. The COSMOS2020 sources which have IRAC channel 3 (5.8 $\mu$m) detections as shown as green triangles. The orange diamonds are all our MIRI-detected sources in this field. Note that we have $>$3$\times$ more MIRI 5.6 $\mu$m detections in the field relative to the prior Spitzer IRAC 5.8 $\mu$m detections.

Figure 11: The redshift distribution of our COSMOS2020 matched MIRI 5.6 $\mu$m sources. These are LePhare redshifts from the classic COSMOS2020 catalog.

4.3 The MIRI 5.6 $\mu$m sources in the COSMOS2020 catalog

As already discussed in Section 3.3, we cross-matched the sources within our two COSMOS fields with the COSMOS2020 catalog (Weaver et al., 2022a, b). The overlap between our MIRI detections, the COSMOS2020 sources and their subset with IRAC channel 3 detections is illustrated in Figure 10 in the case of the COSMOS1 field. Note that with these new data more than half of the COSMOS2020 sources now have 5.6 $\mu$m detections. The number of sources with MIRI 5.6 $\mu$m detection in this field is more than 4$\times$ the sources with previous IRAC channel 3 (5.8 $\mu$m) detections. This is illustrative of the fact that until JWST, it was difficult to link the optical and mid-IR source populations. Within our two COSMOS MIRI fields, there are 690 MIRI sources, of which 582 (84%) are matched to COSMOS2020 catalog sources, using a 1$\arcsec$ matching radius (see Section 3.3).

Figure 11 shows the redshift distribution of the MIRI sources with counterparts in the COSMOS2020 catalog. We show the LePhare redshifts; although the distribution is very similar using the EAZY redshifts. In Figure 11, we applied a quality cut of $0.01<|z_{u68}-z_{l68}|/(1+z_{\mathrm{med}})<0.5$, where $u68$, $l68$, and $\mathrm{med}$ denote the upper and lower 68% confidence intervals and the median of the redshift probability density function. Of the 582 matched sources, 520 passed the quality cut and are presented the figure. This histogram shows peaks at $z\sim 1$ and $z\sim 2$, i.e. we find many cosmic noon galaxies.

Up to $z\sim 2.5$, detection is aided by a negative $k$-correction as the observed 5.6 $\mu$m probes up the SED towards the stellar 1.6 $\mu$m bump, which explains the sharper drop in sources at higher redshift. However, we do see a tail out to $z\sim 5$ which corresponds to early massive galaxies.

Figure 12 shows the stellar mass vs. star-formation rate for our COSMOS2020-matched MIRI sources. We find predominantly star-forming main sequence galaxies with stellar masses in the $10^{8}-10^{10}$ M_⊙ range. We are probing the typical low-mass star-forming galaxy at cosmic noon, well below $M^{*}$ ($\approx 10^{11}$M_⊙) at cosmic noon (Adams et al., 2021).

4.4 The nature of MIRI sources not in COSMOS2020

As discussed in Section 3.4, while the bulk of our sources (84%) have counterparts in the COSMOS2020 catalog, 108 MIRI are un-matched. We examined these by eye to remove sources that have likely photometric issues such as being driven by a single hot pixel or being residuals of the stripes we saw prior to background subtraction (see Figure 2). After this visual inspection, we are left with 31 reliable MIRI detections not present in the COSMOS2020 catalog. This translates to a source density of $\sim$7/sq.arcmin.

Figure 12: The stellar mass vs. star formation rate relation for the same sources presented in Figure 11. The overlaid lines are the main sequence relations at the indicated redshifts based on Koprowski et al. (2024).

To understand why strong MIRI detected sources may not be in the COSMOS2020 catalog consider that the COSMOS2020 source detection is based on a $\chi^{2}$-weighted $izYJHK$ image (Weaver et al., 2022a). Figure 13 shows the 5 $\sigma$ limits of the COSMOS2020 catalog in these six bands, as well as the four IRAC channels. In all cases, we take the limits from the deep parts of COSMOS since these correspond to the locations of our two COSMOS fields. Note that our MIRI F560W depth is comparable to the IRAC channel 1 and 2 depths in this field, while the IRAC channel 3 image is approximately 3 magnitudes shallower (see figure 13). We overlay two potential models that could explain non-inclusion in the COSMOS2020 catalog, but detection in MIRI 5.6 $\mu$m. One is a $z=1$ $M_{star}=10^{9}M\odot$ star-forming galaxy with $A_{V}=3$ (the “Low-mass dusty” model). The other is a $z=4$ $M_{star}=5\times 10^{10}M\odot$ star-forming galaxy with $A_{V}=3$ (the “High mass dusty” model). Here we use the CSP models from Maraston (2005), assume an age of 300Myr for both galaxies and use the Calzetti law (Calzetti et al., 1994) for dust attenuation. The first model is included because while, previously, low mass galaxies were typically not considered very dusty, JWST has recently found evidence of a population of $z<2$ extremely dusty dwarfs (with masses $10^{7}-10^{9}M_{\odot}$ and up to $A_{V}\sim 5$ , Bisigello et al. 2023; see also Pope et al. 2023).

Figure 13: We compare the 5 $\sigma$ magnitude limit of our images (in red) with the $izYJHK$+IRAC 5 $\sigma$ magnitude limits for the COSMOS2020 catalog (in grey), based on Weaver et al. (2022a). We suggest two potential explanations for the MIRI detections that are not in the COSMOS catalog: either a low mass intermediate ($z\sim 1$) redshift galaxy or a galaxy closer to $M_{*}$ at $z\sim 3$ or beyond with high dust content. See Section 4.4 for more details on the models.

In Figure 13, our model spectra are chosen to both lie below the limits for source detection in COSMOS2020 but also to have MIRI 5.6 $\mu$m magnitudes on par with those of the sources shown in Figure 14. Both models satisfy these criteria, but they are clearly distinguishable by the IRAC channel 1 brightness. Therefore in Figure 14 we overlay our MIRI images of these 31 sources with the 2 and 3 $\sigma$ contours from the COSMOS IRAC channel 1 image. We used the IRAC channel 1 mosaic combining all available COSMOS data, as produced by Annunziatella et al. (2023). We find that 3-7 of the sources have emission in IRAC channel 1. The range given is to distinguish clean detections from ones affected by blending in the IRAC image. The three clean IRAC detections in particular (ID 1-16, 1-107, 1-110) are all consistent with the low-mass dusty model. However, the bulk of our 31 sources are undetected in these IRAC images suggesting they are higher-mass and higher redshift dusty sources. Further study of these 31 sources, folding on the new COSMOS-Web data (Casey et al., 2023), as well as all other available data in COSMOS is reserved for a separate paper.

Figure 14: $3.3\times$3.3^′′ cutouts of the 31 sources detected in MIRI 5.6 $\mu$m but are not present in the COSMOS 2020 Classic catalog. The red contours represent the 2 and 3 $\sigma$ levels from the COSMOS IRAC channel 1 image from Annunziatella et al. (2023), with an applied Gaussian blur ($\sigma=0.8$). We label each source with our catalog ID and its MIRI 5.6 $\mu$m flux.

5 Summary & Conclusions

In this paper, we present a study of the number counts and source properties of sources detected in eight MIRI F560W images with a combined area of 18.4 sq.arcmin and a 5 $\sigma$ depth of $\approx$0.1 $\mu$Jy. Two of our eight fields overlap with the COSMOS2020 catalog which also allows for a more detailed look at the source population properties. Below we summarize our key findings:

• The pipeline reduced MIRI F560W images required further background subtraction, given strong stripping artifacts. We also note some residual structure seen in the noise, consistent with tree ring artifacts as noted in $JWST$ documentation. The photometric zero point needed to be adjusted since the header keyword in these early images was not using the true measured MIRI PSF.

• Our study includes 8 widely separated fields, allowing us to explore the effects of cosmic variance. We used this to construct field-to-field variation uncertainty on our combined MIRI 5.6 $\mu$m number counts.

• Our number counts have a more pronounced knee, at $\approx$ 2 $\mu$Jy and are $\approx 2\times$ higher than those computed from the CEERS images (Yang et al., 2023). This difference however is only at $\approx$1 $\sigma$ level given our measured field-to-field variation. The observed counts are consistent with the Cowley et al. (2018) SAM predictions around the knee. These models however overpredict the observed counts below the knee of the counts.

• We find 84% of the MIRI sources in the COSMOS fields have counterparts in the COSMOS2020 catalog. They are predominantly at cosmic noon, and often dusty star-forming galaxies with stellar masses well below $M^{*}$.

• We find 31 very red sources that are securely detected in MIRI, but are not in the COSMOS2020 catalog. This population has a source density of $\approx$7/sq.arcmin. They are consistent with being either intermediate redshift ($z\sim 1$) dusty lower mass galaxies (those with a strong IRAC channel 1 detection, 10–20%) or high mass, high-$z$ ($z\gtrsim 4$) dusty galaxies (those with weak or no IRAC channel 1 detections, 80–90%). These MIRI 5.6 $\mu$m only sources will be explored further in a follow-up paper.

All the JWST data used in this paper can be found in MAST: http://dx.doi.org/10.17909/jsqw-mq02 (catalog 10.17909/jsqw-mq02).