Star-forming brightest cluster galaxies at z ∼ 0 . 4 in KiDS. Further studies of cold gas and stellar properties.

,

Simulations predict that a large fraction of the BCG stellar mass is assembled at high redshifts (z ∼ 3−5), while smaller galaxies are later swallowed by the BCG itself: a process that allows the BCGs to grow in mass and size (De Lucia & Blaizot 2007).BCGs indeed double their stellar mass since z ∼ 1 (Lidman et al. 2012), which is consistent with evolution via dry accretion of satellites (Collins et al. 2009;Stott et al. 2011).However, recent studies have proposed that high star formation and large molecular gas reservoirs are present in a number of intermediate-redshift BCGs (Fogarty et al. 2015(Fogarty et al. , 2019;;Castignani et al. 2020aCastignani et al. , 2022a;;Dunne et al. 2021;O'Sullivan et al. 2021), where star formation is fed by rapid gas deposition, possibly induced by recent mergers or interactions with the companions (Castignani et al. 2022a).These results imply strong environmental quenching mechanisms for distant BCGs.Therefore, the standard scenario of an early formation of BCGs, with a peak of star formation activity followed by passive evolution, needs to be reconsidered, at least for the subpopulation of these distant BCGs with significant ongoing star formation activity (i.e., with a star formation rate, SFR, 20 M yr −1 ) and with large molecular (H 2 ) gas reservoirs (M H 2 10 10 M ).
In recent studies, we have started a large campaign exploiting millimeter facilities with the final goal of understanding the impact of the cluster environments in processing the galactic gas reservoirs, with a particular focus on distant (z > 0.3) BCGs, their mass assembly, and their star formation properties (Castignani et al. 2019(Castignani et al. , 2020a(Castignani et al. ,b,c, 2022a,b),b).The advent of ongoing and forthcoming infrared-to-optical surveys such as the Dark Energy Survey (McClintock et al. 2019;Aguena et al. 2021), the Vera Rubin Observatory Legacy Survey of Space and Time (LSST; Ivezić et al. 2019), and Euclid (Euclid Collaboration 2022) will revolutionize the field in the next years, with thousands of clusters up to z ∼ 2, bringing an unprecedented leverage for constraining the evolution of BCGs in clusters (Rhodes et al. 2017).
These surveys also include the multiwavelength Kilo Degree Survey (KiDS; de Jong et al. 2017;Kuijken et al. 2019).In its final release, this imaging survey will cover about 1350 deg 2 of the sky (de Jong et al. 2017;Kuijken et al. 2019), equally divided between an equatorial (KiDS-N) and a southern field (KiDS-S).KiDS exploits VLT Survey Telescope (VST) ugri optical observations and Visible and Infrared Survey Telescope for Astronomy (VISTA) telescope ZY JHK infrared imaging, and it is able to detect sources down to a limiting AB magnitude r ∼ 25.0 (5σ in 2 aperture).Optical-infrared photometry enables the determination of accurate photometric redshifts up to z ∼ 1 (Benítez 2000;Hildebrandt et al. 2021).These are complemented by spectroscopy, as KiDS is also covered by the Galaxy And Mass Assembly (GAMA; Baldry et al. 2010) and the Sloan Digital Sky Survey (SDSS; York et al. 2000) spectroscopic surveys.This rich data set has enabled the detection of a sample of ∼8000 galaxy clusters up to z ∼ 0.8 in KiDS at a signal-to-noise ratio (S /N) > 3 (Maturi et al. 2019;Bellagamba et al. 2019) with the Adaptive Matched Identifier of Clustered Objects (AMICO) cluster finder (Bellagamba et al. 2018), which searches for overdensities of galaxies using photometric redshifts with a matched filtering.
In this work, we study three star-forming BCGs that we selected as those among the most star-forming BCGs in KiDS.We present our analysis resulting from new observations of the three BCGs in the first three CO transitions with the Institut de Radioastronomie Millimétrique (IRAM) 30 m telescope.We follow a similar approach to that of our recent study, Castignani et al. (2022a), hereafter denoted C22a, where we investigated the molecular gas and star formation properties of three additional BCGs drawn from the same parent sample.In C22a, we were able to double the number of distant BCGs with clear detections in at least two CO lines based on new molecular gas observations.To explain the global stellar, gas, and star formation properties of the BCGs, we suggested that a substantial amount of the molecular gas might have been accreted by the BCGs, possibly via accretion by the nearby cluster core companions, but might still not have been efficiently converted into stars.Alternatively, the observed molecular gas reservoirs could originate from the cooling of the hot X-ray atmospheres around the BCGs.We also showed, for the first time for distant BCGs, that gas excitation is correlated with the specific star formation rate (sSFR) and the star formation efficiency at which molecular gas is consumed, as was previously found only in other types of star-forming galaxies.Highly excited gas is present only in cool-core and highly star-forming (SFR 100 M yr −1 ) BCGs, where the low-entropy cool cores of clusters sustain high levels of star formation and favor the condensation and excitation of cold gas in the BCGs (Castignani et al. 2020a).Together, these results suggest that environmentally driven mechanisms play an important role in processing the molecular gas that feeds star formation in distant BCGs.In order to better understand the physical mechanisms that cause this processing, it is thus essential to increase the still limited sample of distant BCGs with detections at different molecular gas transitions, which probe different gas densities, and possibly different degrees of environmental processing.
The paper is structured as follows.In Sect. 2 we introduce the BCGs that are the subject of this work, in Sect. 3 we describe the molecular gas observations and analysis, in Sect. 4 we present the results, and in Sect. 5 we summarize them and draw our conclusions.Stellar mass and SFR estimates reported in this work for the three targeted BCGs rely on the Chabrier (2003) initial mass function (IMF).Magnitudes are reported in the AB system.Throughout this work, we adopt a flat Λ-cold dark matter (ΛCDM) cosmology with matter density Ω m = 0.30, dark energy density Ω Λ = 0.70, and Hubble constant h = H 0 /(100 km s −1 Mpc −1 ) = 0.70.

Three star-forming
BCGs at z ∼ 0.4

Selection of the BCGs in KiDS
One of the main goals of this work is to search for molecular gas reservoirs in distant BCGs and to characterize their stellar and star formation properties.To do this, similarly as C22a, we started by considering a parent sample of 1484 BCGs with spectroscopic redshifts in the equatorial KiDS-N area.Of these, 684 BCGs are at z > 0.3.These BCGs were selected by Radovich et al. (2020) using the Data Release (DR) 3 of KiDS (de Jong et al. 2017).As we are interested in actively star-forming and distant BCGs, we further limited ourselves to the BCGs with clear counterparts in the all-sky data release of the Wide-field Infrared Survey Explorer (WISE; Wright et al. 2010) at 22 µm in the observer frame (i.e., the W4 channel).We found only six spectroscopically confirmed BCGs at z > 0.3 with WISE W4 S /N > 3.3.The three of these with the highest values, S /N ∼ 4.1−6.5, were analyzed in C22a, while in this study, we consider the remaining three, which have WISE W4 S /N ∼ 3.4−3.8.The chosen redshift threshold of z > 0.3 allows us to focus on distant BCGs and on their evolution, while the cold gas content of more local BCGs, including those of Abell clusters, has been previously and extensively studied (e.g., Edge 2001;Salomé & Combes 2003;Olivares et al. 2019;Rose et al. 2023).

Stellar and star formation properties
In the following, we investigate the multiwavelength properties of the three BCGs of this work.In Fig. 1 we show the colorcomposite images of the three galaxies.KiDS 1433 is bulge dominated, while KiDS 0837 and KiDS 1415 appear to be have a more disturbed morphology and clumpy substructures.The visual inspection of the images suggests that our targets are in dense environments.Similarly to the KiDS BCGs studied in C22a, the three BCGs have at least one companion within a projected separation of ∼4.2 arcsec (i.e., ∼23 kpc at z = 0.4).

Spectral energy distributions
In Fig. 2 we show the ultraviolet to far-infrared spectral energy distributions (SEDs) of the three sources, taken from the GAMA database.The SED modeling reported in the panels of the figure corresponds to the MAGPHYS (da Cunha et al. 2008) fits by Driver et al. (2018), which are publicly available in GAMA DR3.Photometric data include GALEX (Martin et al. 2005;Morrissey et al. 2007) in the ultraviolet, SDSS (York et al. 2000) in the optical, and the VISTA Kilo-degree Infrared Galaxy Survey (VIKING; Edge et al. 2013), WISE (Wright et al. 2010), and Herschel-ATLAS (Eales et al. 2010;Valiante et al. 2016) in the near-to far-infrared.
In Table 1 we summarize the main properties of the BCGs, including their stellar, star formation, and dust properties derived from the SED fits (Driver et al. 2018).The stellar masses of all three BCGs exceed 10 11 M , which confirms that they are very massive galaxies.Furthermore, all three BCGs have far-infrared emission from Herschel, well modeled by a dust component.This supports the scenario that the observed infrared emission is mainly due to star formation, with minimum AGN contamination.Dust masses and luminosities are in the range log(M dust /M ) 8.6−9.4 and log(L dust /L ) 11.4−11.6,respectively, which is typical of luminous infrared galaxies (LIRGs).The SFR estimates based on the ultraviolet to far-infrared SED modeling are high and lie in the range ∼(15−24) M yr −1 .These are consistent with but higher than the SFRs of main-sequence (MS) galaxies of similar masses and redshifts (Speagle et al. 2014).The only exception is KiDS 0837, which has the lowest SFR, a factor of 2 lower than the average SFR of MS galaxies at its redshift.
Interestingly, the stellar masses, dust luminosities, and SFRs of the three BCGs of our sample are similar to those of the three star-forming z ∼ 0.4 KiDS BCGs we studied recently in C22a.The six sources were selected from the KiDS DR3 sample of 684 BCGs with spectroscopic redshifts z > 0.3 (Sect.2.1).They represent 0.8% of this parent sample and are among the most star-forming BCGs in KiDS.As already pointed out in C22a, this rare population of BCGs experiences stellar mass assembly and is thus caught in a special phase of their evolution (see, e.g., Castignani et al. 2020a,b, for similar BCGs).These BCGs are thus the intermediate-z counterparts of local star-forming BCGs (>40 M yr −1 ) such as the famous Perseus A and Cygnus A (Fraser-McKelvie et al. 2014).

Line diagnostics
We now investigate the star formation properties of the BCGs using the spectra found in the GAMA database.We performed a similar analysis to that presented in C22a, as outlined in the following.The SFR is directly proportional to both Hα line and [O II] forbidden-line doublet (3726 Å, 3729 Å) luminosities L Hα and L [O II] , respectively (Kennicutt 1998).We therefore used these relations to estimate the SFR, calibrated to a Chabrier (2003) IMF.We also corrected for dust attenuation using the empirical relation by Garn & Best (2010), which relates dust attenuation (A Hα ) to stellar mass as follows: where M = log(M /M ) − 10.With these conventions, the Hαand [O II]-based SFRs can be expressed as follows (see, e.g., Gilbank et al. 2010;Zeimann et al. 2013;Old et al. 2020):  where r lines is the ratio of extinguished [O II] to Hα fluxes.We found GAMA spectra for all three BCGs considered in this work, taken with the AAOmega/2dF multifiber spectrograph on the 3.9 m Anglo-Australian Telescope.Gaussian fits to several lines are also available (Gordon et al. 2017).Unfortunately, Hα is redshifted to wavelengths longer than 8890 Å for all three sources, and they are therefore beyond the wavelength range of 3750 Å λ 8850 Å that is covered by the spectrograph (Hopkins et al. 2013).
However, the [O II] doublet is detected in the GAMA spectra of all three BCGs.Following Gilbank et al. (2010) With these fluxes, Eq. ( 3) yields an SFR [O II] = (5.5 ± 1.7), (21.9 ± 7.4), and (25.1 ± 7.6) M yr −1 for KiDS 0837, 1415 and 1433, respectively, where we assumed r lines = 0.5, following Gilbank et al. (2010).With the only possible exception of KiDS 0837, which has the lowest SFR, we find a good agreement within the reported uncertainties when we compare the [O II]based SFRs with those obtained by SED fits (Table 1), which are used throughout this work.This is remarkable considering that discrepancies, even up to a factor of ∼3, are common when different estimates of the SFR are compared (e.g., Calzetti 2013).

WISE infrared colors
Finally, we investigated the star formation properties and the AGN contamination in the infrared by inspecting the position of our sources in the infrared W1 − W2 vs. W2 − W3 color-color WISE diagram (Jarrett et al. 2017).In this diagram, sources with W1 − W2 0.8 mag can be considered to be AGN, while lower values of W1 − W2 are associated with other types of sources.The W1 − W2 colors of the three KiDS BCGs of this work lie in the range between (0.36−0.55) mag and are therefore not classified as AGN.On the basis of their WISE colors, both KiDS 1415 and 1433 are indeed classified as starbursts, while KiDS 0837 falls in the region of galaxies with SFRs below the starburst level.However, the last BCG has only a 3σ lower limit of W3 > 11.2 mag.This implies W2 − W3 < 2.7, so that the spheroid class is still possible for the source.Overall, these findings strengthen the selection of our three targets as star-forming galaxies and suggest that any possible contamination from a circum-nuclear dusty torus at infrared wavelengths is likely negligible, which is also supported by the generally good agreement between the SED-and line-based SFRs.Interestingly, C22a found similar line-based SFRs and infrared colors A139, page 4 of 15   3) and ( 4) CO(J → J − 1) transition and observer frame frequency, ( 5) CO(J → J − 1) velocity-integrated flux, (6) S/N of the CO(J → J − 1) detection, ( 7) FWHM of the CO(J → J − 1) line; ( 8) CO(J → J − 1) velocity-integrated luminosity, (9) redshift derived from the CO(J → J − 1) line, and (10) on-source integration time (double polar).Absent values are denoted with the long dash.For KiDS 0837 and KiDS 1415, the reported upper limits are at 3σ, estimated at a resolution of 300 km s −1 .
for the additional three star-forming intermediate-z KiDS BCGs that come from our selection.These results thus show that the six BCGs together belong to a homogeneous population of BCGs in terms of stellar, star formation, and dust properties.However, as we discuss in Sect.4.3, their molecular properties are more heterogeneous.

Millimeter IRAM 30 m observations
We observed the three KiDS BCGs using the IRAM 30 m telescope at Pico Veleta in Spain.The observations were carried out on December 15−21, 2021 (ID: 176-21; PI: G. Castignani).We used the Eight Mixer Receiver (EMIR) to observe CO(J → J−1) emission lines from the BCGs, where J is a positive integer denoting the total angular momentum.For each galaxy, the specific CO(J → J − 1) transitions were chosen to maximize the detection likelihood in terms of the ratio of the predicted signal to the expected root mean square (rms) noise.The E090, E150, and E230 receivers, operating at ∼1−3 mm, offer 4 × 4 GHz instantaneous bandwidth covered by the correlators.For the KiDS 0837 and KiDS 1415 BCGs, we used the two E090 and E150 receivers simultaneously to observe CO(1→0) and CO(2→1), redshifted to 3 mm and 2 mm in the observer frame, respectively.For all three targets, we also used the E090 and E230 receivers to simultaneously target the CO(1→0) and CO(3→2) lines, redshifted to 3 mm and 1 mm, respectively.We note that we did not observe KiDS 1433 in CO(2→1) as the line is redshifted to 170.189 GHz in the observer frame, which is too close to the 183 GHz water vapor line.The line is thus prohibited because of low atmospheric transmission.
The IRAM 30 m half-power beam width (HPBW) is ∼16 arcsec λ obs 2 mm (Kramer et al. 2013), where λ obs is the observer frame wavelength.All targets are thus unresolved by our observations, assuming a CO-to-optical size ratio ∼0.5 (Young et al. 1995); see Fig. 1.The wobbler-switching mode was used for all the observations to minimize the impact of atmospheric variability.The Wideband Line Multiple Autocorrelator (WILMA) was used to cover the 4 × 4 GHz bandwidth in each linear polarization.We also simultaneously recorded the data with the fast Fourier Transform Spectrometers (FTS) as a backup at 200 kHz resolution.
We encountered variable weather conditions during our observations.Overall, we had average system temperatures typical of the winter season and equal to T sys 107−116 K, 228−278 K, and 183−224 K at 3, 2, and 1 mm, respectively.As calibrators, we used the quasars PKS 1253−055 and OJ 287 for pointing, and PKS 1253−055 for the focus.We summarize our observations in Table 2.The data reduction and analysis were performed using the class software of the gildas package2 .
Our results are presented in Sect.4, and in the following, we describe our analysis of the IRAM-30 m spectra for the three targeted BCGs.

Molecular gas properties
For KiDS 1433, our observations yield clear detections of CO in the first and third transitions at an S /N 5.0−5.5.The corresponding spectra are shown in Fig. 3.We fit the total emission using a single Gaussian for each of the CO(1→0) and CO(3→2) spectra of KiDS 1433.As the CO(1→0) and CO(3→2) lines of KiDS 1433 are double-peaked, we additionally fit the total emission for each spectrum using a combination of two Gaussian curves that are overlaid on the spectra in Fig. 3.In Table 2 we report the best-fit results for the single-Gaussian fits (denoted Total) and for the double-Gaussian fits (labeled v > 0 and v < 0).For the last ones, v > 0 (v < 0) denotes the fit to the CO emission that peaks at positive (negative) velocity with respect to the BCG redshift.As outlined in Table 1, we observe an overall agreement between the single-and double-Gaussian fits in terms of velocity-integrated flux, full width at half maximum (FWHM), and S/N.Our CO-based redshift measurement also agrees well with the BCG redshift that is estimated with optical spectral lines.The observed double-peaked CO emission is probably likely due to the KiDS 1433 BCG alone, and we further discuss the underlying structural parameters in Sect.4.7.However, the BCG has two nearby companions at a projected separation of ∼4.2−4.7 arcsec (see Fig. 1).They have high AMICO cluster membership probabilities of 0.93−0.96and are therefore photometrically selected cluster members.These BCG companions may contribute to the CO(1→0) and CO(3→2) emission, which may accordingly result in the observed double-horn profile.Nevertheless, this appears to be unlikely as the efficiency is highly suppressed for CO(3→2), down to 54% at most, given the angular separations of the BCG companions.For CO(1→0), the efficiency suppression is only marginal, at a level of 93%.However, it is worth noting that the positive and negative peaks have similar height, and each of them similarly contributes to (46±11)% of the total emission for both CO(1→0) and CO(3→2).Similarly, the v > 0 peaks and the v < 0 peaks have similar velocities in the two panels of Fig. 3.These results indicate a scenario in which the double-peak emission in the CO(1→0) and CO(3→2) spectra has the same physical origin, and we interpret it as CO emission from the BCG KiDS 1433 itself.However, higher angular resolution observations, for example, with ALMA or NOEMA interferometers, are needed to spatially resolve the molecular gas reservoir and firmly test the proposed scenario.
For the other two targeted BCGs, KiDS 0837 and 1415, we did not detect any CO emission in their spectra.We therefore assigned upper limits to the CO line fluxes as follows.We first removed the baseline in each spectrum with a linear fit.We then estimated the rms noise level for the antenna temperature (Ta * ).Using the standard conversions outlined below, we converted the rms into a 3σ upper limit to the integrated CO flux at a resolution of 300 km s −1 , which is the typical FWHM for massive galaxies and BCGs in particular (Edge 2001;Dunne et al. 2021).
The results of our analysis are reported in Table 2 for all three BCGs we targeted in this work.We applied standard efficiency corrections to convert (i) the antenna temperature Ta * into the main beam temperature T mb , and then (ii) T mb into the corresponding CO line flux, with a conversion factor of 5 Jy K −1 .In particular, we adopted the following efficiency corrections: T mb /T a * = 1.2, 1.4 and 2.0 for ∼3, 2, and 1 mm observations, respectively 3 .
3 https://www.iram.es/IRAMES/mainWiki/Iram30mEfficiencies We then derived the CO(J → J − 1) luminosity L CO(J→J−1) , in units of K km s −1 pc 2 , from the velocity-integrated CO(J → J − 1) flux S CO(J→J−1) ∆v, in units of Jy km s −1 .To do this, we used Eq. ( 3) from Solomon & Vanden Bout (2005), where ν obs is the observer frame frequency in GHz of the CO(J → J − 1) transition, D L is the luminosity distance in Mpc, and z is the redshift of the BCG.
We then converted the velocity-integrated CO(1→0) line fluxes into total H 2 gas masses, or into their upper limits.The CO(1→0) transition is preferred to higher-J transitions as it does not require any assumption on the excitation ratio r J1 = L CO(J→J−1) /L CO(1→0) .
The three KiDS BCGs have similar sSFRs within the range sSFR (0.5−2.2) × sSFR MS (see Table 1).They are all lower than the value of 3 × sSFR MS , below which the galaxy sSFR is within the characteristic scatter of the MS.To estimate H 2 gas masses, we therefore assumed a Galactic CO-to-H 2 conversion factor of α CO = 4.36 M (K km s −1 pc 2 ) −1 , which is commonly used for star-forming galaxies at the MS.The use of a single α CO conversion factor also allowed us to perform a homogeneous comparison in terms of gas content (see Sect. 4).
Last, we used the H 2 gas mass estimates of the BCGs to estimate a number of related quantities, or upper limits.The first quantity is the depletion timescale τ dep = M H 2 /SFR, which is the characteristic time in which the molecular gas reservoirs are depleted.To compute τ dep , we used the SED-based SFR estimates reported in Table 1.Similarly, we also estimated the molecular-gas-to-stellar-mass ratio, M H 2 /M and the ratio of the molecular gas to the dust masses.In addition, as a comparison, we computed the depletion time τ dep,MS and the molecular gasto-stellar mass ratio M H 2 M MS for MS galaxies in the field with a redshift and stellar mass equal to those of our BCGs, by using the empirical prescriptions by Tacconi et al. (2018), calibrated to the CO-to-H 2 conversion factor of α CO = 4.36 M (K km s −1 pc 2 ) −1 we used here.
The resulting molecular gas properties for the three KiDS BCGs are summarized in Table 3.As KiDS 1433 BCG is detected in both CO(1→0) and CO(3→2), we also report the associated r 31 excitation ratio in the table.
For KiDS 0837 and KiDS 1415, the reported upper limits are at 3σ.

Results
One of the main results of this work is that with the new CO(1→0) and CO(3→2) detections of KiDS 1433, we notably increase the still limited sample of distant BCGs with clear detections in two distinct CO lines.To the best of our knowledge, only five additional BCGs are detected in multiple CO transitions, all of them at intermediate redshifts z ∼ 0.4.They are the star-forming BCGs KiDS 0920, 1220, and 1444 (C22a), and the BCGs RX1532 (Castignani et al. 2020a) and M1932 Fogarty et al. (2019) from the CLASH survey (Postman et al. 2012).In the following, we present and further discuss our results in detail.

Comparison sample
We compared the molecular gas properties of our targeted BCGs with those of existing samples of cluster galaxies with molecular gas observations.To do this, we considered the sample of z < 2 cluster galaxies with M > 10 10 M described in Sect.4.1 of C22a and references therein, which is an update of the compilation reported in Tables A.1 and A.2 of Castignani et al. (2020b), to which we refer for further details.The only difference with respect to the comparison sample used in C22a is as follows.In C22a, we mostly considered cluster galaxies with detections in CO.However, for two out of the three BCGs of this work we are only able to set upper limits to their molecular gas mass.Therefore, we now considered the full sample of Dunne et al. (2021), consisting of distant BCGs within the SpARCS survey.
More specifically, we included the BCGs in their sample in our analysis with CO detections and 3σ upper limits.The full list of comparison galaxies we considered comprises 97 cluster galaxies over 55 clusters at z ∼ 0.2−2 and stellar masses in the range log(M /M ) 10−12.They include 40 BCGs that we briefly outline in the following.
We included the gas-rich Cl J1449+0856 BCG at z = 1.99 (Strazzullo et al. 2018), with M H 2 = 7.5 × 10 10 M (Coogan et al. 2018), and the SpARCS 104922.6+564032.5BCG at z = 1.7 (Webb et al. 2017;Barfety et al. 2022), with an upper limit of M H 2 < 1.1 × 10 10 M to the molecular gas mass.The recent molecular gas survey by Dunne et al. (2021) with ALMA adds another 30 SpARCS BCGs: seventeen of them are detected in CO(2→1) and S /N > 3, and we report for the remaining 13 3σ upper limits to their H 2 gas mass here.These BCGs span the redshift range z ∼ 0.2−1.2 and their molecular gas masses lie in the range M H 2 ∼ (10 10 −10 11 ) M , or upper limits in the range M H 2 (0.3−4.8) × 10 10 M .
Last, we included the three BCGs KiDS 0920, 1220, and 1444 at z ∼ 0.4, with gas masses in the range M H 2 ∼ (5−15) × 10 10 M , as inferred from CO(1→0) detections (C22a).In Sect.4.2 we use this sample of distant cluster galaxies observed in CO, including BCGs, to compare stellar mass, SFR, molecular gas content, and depletion time with the three BCGs of this work.
In addition to the ∼100 cluster galaxies with observations in CO listed above, we considered a second compilation of field galaxies with molecular gas observations in the first and third CO transitions that we use in Sect.4.5 to investigate the excitation ratios.Analogously to what we discussed in Sect.4.3 of C22a, this comparison sample of sources with observations in multiple CO transitions includes 5 distant (2.5 < z < 2.7) galaxies that are part of the ALMA Spectroscopic Survey in the Hubble Ultra Deep Field (ASPECS; Boogaard et al. 2020), and 61 local galaxies and AGN from Lamperti et al. (2020), which are drawn from the xCOLD GASS survey (Saintonge et al. 2017) and the BAT AGN Spectroscopic Survey (BASS)4 .

Gas fraction, star formation, and depletion time
In this section, we place the molecular gas, stellar, and star formation properties of our targeted BCGs in a broader context.To do this, we used the compilation of ∼100 cluster galaxies with M and SFR estimates as well as observations of their molecular gas, as outlined in Sect.4.1.
Figure 4 (left panel) shows the ratio of molecular gas to stellar mass M H 2 /M as a function of redshift.In the right panel of Fig. 4, we instead show the M H 2 /M ratio as a function of the specific SFR, both normalized to the MS. Figure 5 displays the SFR, M H 2 /M , and the depletion time (τ dep ), all normalized to their MS values, as a function of the stellar mass.
In both figures, we report data of distant cluster galaxies with molecular gas observations.They include the three targeted Fig. 4. Molecular gas properties of distant BCGs and cluster galaxies observed in CO.Left: evolution of the molecular gas-to-stellar mass ratio as a function of the redshift for cluster galaxies at 0.2 z 2 observed in CO.The solid green curve is the scaling relation found by Tacconi et al. (2018) for field galaxies in the MS and with log(M /M ) = 11, which corresponds to the median stellar mass of all sources in the plot.The dashed green lines show the statistical 1σ uncertainties in the model.Right: molecular gas-to-stellar mass ratio plotted against the sSFR for the cluster galaxies, normalized to the corresponding MS values using the relations for the ratio and the SFR by Tacconi et al. (2018) and Speagle et al. (2014), respectively.The dot-dashed green lines show different depletion times in units of the depletion time at the MS.The points in both panels show cluster galaxies with molecular gas observations.The color code is reported at the top left of each panel.KiDS BCGs of this work are reported as filled red circles, and open circles refer to KiDS BCGs in C22a.Additional distant BCGs from the literature are also highlighted.We refer to the text for further details.
BCGs of this work (filled red circles), the three KiDS BCGs of our recent C22a work (open red circles), and the compilation of ∼100 cluster galaxies of our comparison sample, described in Sect.4.1.In the latter, BCGs from the literature are highlighted, and the remaining cluster galaxies in our comparison sample are shown in background (gray dots).To allow a homogeneous comparison, we used α CO = 4.36 M (K km s −1 pc 2 ) −1 for all data points and the overplotted Tacconi et al. (2018) relations in both figures.
Our targeted BCGs have high SED-based stellar masses in the range of log(M /M ) ∼ (11.0−11.6);see Table 1.As illustrated in Fig. 5, these BCGs fall in the upper tail of the M distribution for cluster galaxies, which is populated by BCGs in our comparison sample.These results place the targeted sources of this work among the most massive cluster galaxies at intermediate redshifts, strengthening their selection as BCGs (Radovich et al. 2020), which we further discuss in Sect.4.6.Furthermore, Figs. 4 and 5 both show that KiDS 1433 is very rich in molecular gas having M H 2 6 × 10 10 M and a H 2 -to-stellar mass ratio of M H 2 /M 0.32 (see also Table 3).These values agree with those we found in C22a for similarly star-forming z ∼ 0.4 BCGs KiDS 0920, 1220, and 1444, suggesting that the four BCGs together belong to the same population of star-forming and gas-rich BCGs.Interestingly, the other two z ∼ 0.4 BCGs that we targeted in this work, that is, KiDS 0837 and 1415, have similar star formation levels as the other four KiDS BCGs, well within the MS; see Fig 4 (right panel) and Fig. 5 (top panel).However, KiDS 0837 and 1415 BCGs were not detected in CO with our observations, and we place upper limits of M H 2 /M 0.07 and 0.23, respectively, which are consistent with the most stringent upper limits found for distant BCGs in previous studies (e.g., Castignani et al. 2020a,b,c;Dunne et al. 2021).For a direct comparison, we also refer to the upper limits from the literature that are reported in Fig. 4.Under the assumption that the total infrared emission is a tracer of ongoing star formation fueled by the molecular gas reservoir (e.g., Solomon et al. 1997;Tacconi et al. 2018), a possible interpretation is that our millimeter observations are just below the CO detection limit for KiDS 0837 and 1415 BCGs.Their total infrared fluxes are 2.5 and 1.8 × 10 12 erg cm −2 , respectively, as inferred from the total SED-based dust luminosities reported in Table 1.Conversely, the remaining four KiDS BCGs detected in CO all have slightly higher infrared fluxes, in the range of ∼(2.6−4.7)× 10 12 erg cm −2 .However, based on our CO upper limits, even with deeper observations (e.g., with NOEMA or with the ALMA interferometers), we expect KiDS 0837 and 1415 to have an H 2 amount at the level of the MS at most, which is significantly lower than that we infer for the four KiDS BCGs with CO detections, whose H 2 gas masses lie well above the MS levels (Fig. 5, central panel).

Star-forming BCGs in KiDS. Some thoughts on gas-rich versus gas-poor sources
To summarize, the six KiDS BCGs considered in this and our previous work (C22a) have similar levels of star formation activity within the MS.However, as illustrated in Fig. 5, those with CO detections tend to have larger H 2 gas reservoirs and longer depletion times (τ dep ), normalized to the MS values, than those found in the KiDS BCGs with only upper limits in CO.
A139, page 8 of 15 Similarly, there is also a dichotomy in terms of the H 2 -to-dustmass ratio.While the four KiDS BCGs with CO detections have high ratios M H 2 /M dust (170−300), the two BCGs with only upper limits in CO have much lower ratios, well below the value M H 2 /M dust 100, which is typical of distant star-forming galaxies (Scoville et al. 2014(Scoville et al. , 2016;;Berta et al. 2016).
To explain the high values of τ dep , M H 2 /M , and M H 2 /M dust observed in the KiDS BCGs with CO detections, including KiDS 1433 that is the subject of this work, we suggest that they are observed in a peculiar and rare phase of the BCG evolution, in which a substantial amount of the H 2 gas has recently been accreted, but still not efficiently converted into stars.As has been suggested in C22a, the accretion may have occurred via interaction with the cluster companions.The other star-forming BCGs (KiDS 0837 and 1415), with only upper limits in CO, may be at a different and possibly later stage of the BCG evolution, where the H 2 gas reservoir has been exhausted or is undergoing rapid consumption at a timescale τ dep (1−2) Gyr, while star formation is still ongoing.These results suggest that molecular gas depletion precedes star formation quenching in intermediate-z and star-forming BCGs.This result agrees with what Castignani et al. (2022c) found locally for early-type galaxies in the field of the Virgo cluster.
Furthermore, KiDS 0837 and KiDS 1415 BCGs appear to have clumpy substructures (Fig. 1), which is partially at odds with the other KiDS BCGs with detections in CO (see also Fig. 1 of C22a).Therefore, it might be that compact components and clumpy morphology favor the rapid exhaustion of molecular gas and ultimately help to quench the BCGs.This scenario agrees with that proposed in Castignani et al. (2020b), who studied a sample of distant star-forming BCGs mostly with only upper limits in CO and hence short depletion times.By combining molecular gas observations with a Hubble Space Telescope morphological analysis, we speculated that compact components and clumpy morphology may favor the rapid exhaustion of molecular gas and ultimately help to quench the BCGs.This scenario might also be applicable to the KiDS 0837 and KiDS 1415 BCGs.

CO-to-H 2 conversion and metallicities
In the previous section, we discussed an environment-driven scenario that may explain the large molecular gas reservoirs that are observed in some distant star-forming KiDS BCGs.However, as already pointed out in C22a, a lower α CO conversion factor than the Galactic one, as is usually used for starbursts, would lead to lower M H 2 and τ dep .
Therefore, we reconsidered the H 2 gas masses of all six KiDS BCGs studied in this work and in C22a.To do this, we estimated the CO-to-H 2 conversion factor using existing metallicitydependent prescriptions.In particular, we adopted estimates of the metallicity taken from GAMA DR3, which come as outputs of the SED fits with MAGPHYS by Driver et al. (2018, see Fig. 2).We then used the metallicity values to estimate the CO-to-H 2 conversion factor using Eq. ( 4) of Tacconi et al. (2018), which corresponds to the geometric mean of the α CO prescriptions reported in Bolatto et al. (2013) and Genzel et al. (2012Genzel et al. ( , 2015)).Similarly to Tacconi et al. (2018), we relied on the relation log Z = 12 + log(O/H), which is the metallicity on the Pettini & Pagel (2004) scale, and the solar metallicity was assumed to be equal to log Z = 8.67 (see, e.g., Asplund et al. 2009, for a review).
value and the rms dispersion around the median.These results imply that our choice of adopting a Galactic α CO is reasonable for our targeted BCGs, and thus the M H 2 estimates do not appear to be biased.

Excitation ratios
Motivated by the results presented above, we now further investigate the differences of the targeted BCGs with respect to those studied in previous works by means of the excitation ratio.To the best of our knowledge, there are only a few distant BCGs with CO detections in different transitions, all at intermediate redshifts.These are RX1532 (z = 0.361; Castignani et al. 2020a) and M1932 (z = 0.353; Fogarty et al. 2019) from the CLASH survey, with high levels of star formation (SFR 100 M yr −1 ).
In addition, there are other four BCGs from KiDS at redshifts between z = 0.32−0.44 and with lower levels of star formation (SFR 22−34 M yr −1 ).These are the KiDS 0920, 1220, 1444, and 1433 BCGs.The first three were the subject of our recent C22a study, and for the present analysis, we included KiDS 1433, which we targeted in this work.As the sources are all detected in CO(1→0) as well as in CO(2→1), CO(3→2), or CO(4→3), we considered the r 21 , r 32 , and r 43 excitation ratios whenever available.
Similarly to the analysis presented in C22a, in Fig. 6 we show the excitation ratios r J1 of the BCGs against the sSFR (left) and the depletion time (right).For the latter, we used α CO = 4.36 M (K km s −1 pc 2 ) −1 for all BCGs to allow a homogeneous comparison.KiDS 1433 of this work has similar properties in terms of excitation ratio, sSFR, and depletion time as those of the similarly star-forming and gas-rich KiDS BCGs considered in C22a, whose excitation ratios lie in the range r 31 ∼ 0.1 ± 0.3.On the other hand, the CLASH BCGs (M1932 and RX1532), which have higher levels of star formation activity than the KiDS BCGs, also have more excited gas reservoirs, with higher values of r 31 ∼ 0.8−0.9, as well as lower depletion times and higher sSFR values.
As all these six BCGs, including KiDS 1433 from this work, have observations in both CO(1→0) and CO(3→2), we include in Fig. 6 the ASPECS, BASS, and xCOLD GASS galaxies introduced in Sect.4.1 for a comparison, similarly to what we did in C22a.When we consider the r 31 values of the BCGs (KiDS and CLASH) and the comparison galaxies together, we find that the excitation ratio is correlated with sSFR and anticor-related with the depletion time at a level of 3.2σ (p-value = 1.62 × 10 −3 ) and 4.8σ (p-value = 1.26 × 10 −6 ), as found with a Spearman test, respectively.A correlation at the level of 2.8σ (p-value = 4.8 × 10 −3 ) is also found when considering the KiDS and CLASH BCGs alone in the r 31 vs. sSFR plane.These significances are similar as but also slightly higher than those found by a similar analysis in C22a, where KiDS 1433 was not included.
These results strengthen the scenario proposed in our previous studies (Castignani et al. 2020a(Castignani et al. , 2022a) such that while similarly large amounts of molecular gas reservoirs, with M H 2 (0.5−1.4) × 10 11 M and M H 2 /M 0.1−0.6, are often found in star-forming (SFR 20 M yr −1 ) BCGs at intermediate-z, only the most star-forming (SFR 100 M yr −1 ) and cool-core BCGs such as RX1532 and M1932 are able to excite the gas at the highest levels, for which the cool-core cluster environments likely favor the condensation and the inflow of gas onto the BCGs themselves (Castignani et al. 2020a).Low excitation values r J1 0.4 are instead found for the BCGs with high τ dep 2 Gyr such as KiDS 1220 and 1433, which have double-horn CO spectra that are indicative of a low gas concentration (see also Sect.4.7).These results are reasonable, as in the regime of low star formation efficiency, that is, high τ dep , a large amount of gas might be located far out in the disk and not involved in the star formation, as is the case for the nearby giant elliptical Centaurus A (Salomé et al. 2017), for example.
Alternatively, we recall that the IRAM 30 m HPBW is equal to ∼24 and 10 arcsec for the CO(1→0) and CO(3→2) transitions, respectively.Similarly to the three KiDS BCGs of C22a, KiDS 1433 of this work has nearby companions within 10 arcsec in projection.These are clearly visible in Fig. 1 and might contribute to the observed CO(1→0) emission, while the contamination is likely negligible for CO(3→2), as we further discussed in Sect.3.2.The relatively low r 31 values of the KiDS BCGs might therefore be explained when the BCG companions are gas rich.High-resolution interferometric observations are needed for a conclusive answer.

Color-magnitude plots
In this section, we investigate the cluster galaxy population of the three KiDS clusters of this work by means of color-magnitude plots, with a particular focus on the BCGs.We performed an analysis similar to that presented in C22a for the additional three KiDS BCGs with molecular gas observations.We started by selecting sources in the KiDS DR3 within a separation of 2 arcmin from each of the three targeted BCGs.
In Fig. 7 we report the corresponding g − r versus r plots, where g and r are AB apparent magnitudes, obtained with the OmegaCAM camera on the VST (see Sect. 2.1).At the redshifts of our three targeted BCGs, the rest frame 4000 Å break is redshifted to ∼(5500−5600) Å between the g and r filters, whose effective wavelengths are 4800 Å and 6250 Å, respectively.The chosen filters are thus optimal to detect the red sequence in the considered clusters.In the color-magnitude plots, the BCGs are denoted with green stars.We also distinguish field galaxies (gray points) from photometric cluster members (blue points), which were selected as galaxies with cluster membership probabilities greater than 50%, as was further discussed in Sect 5.4 of C22a.Using both GAMA and SDSS spectroscopic databases, we did not find any other spectroscopic cluster member in addition to the BCG themselves.
Visual inspection of the color-magnitude plots shows that our KiDS BCGs are the brightest of the selected cluster members.The one exception is a photometrically selected member with A139, page 10 of 15 For the KiDS 0837 cluster, we detect a strong clustering of red sources at g − r ∼ 1.5.This suggests that the KiDS 0837 cluster core is sufficiently mature and its galaxy population is well evolved in terms of color properties.For the other two clusters, the region of the plot around g − r ∼ 1.5 is instead only weakly populated by cluster members.However, cluster members with bluer colors g − r ∼ 1.0−1.3,similar to those of the associated BCGs appear to lie there as well, which may indicate that ongoing star formation activity is not limited to the BCGs.The BCGs KiDS 0837, 1415, and 1433 are systematically blue, that is, with relatively low values of g − r = 1.4,1.0, and 1.2, respectively.These values are indeed ∼(0.3−0.5)mag lower than the median value of g − r = 1.70 +0.11   −0.255 , which is associated with 0.3 < z < 0.5 BCGs in KiDS DR3 (Radovich et al. 2020).These results generally agree with the fact that we selected the BCGs to be star-forming so that they are caught in a rare phase of their evolution.
These results encouraged us to perform a red-sequence modeling, similarly to previous studies (Kotyla et al. 2016;Castignani et al. 2019).We followed the same approach as was adopted in Sect.5.3 of C22a, to which we refer for further details.In particular, we used the tool Galaxy Evolutionary Synthesis Models (GALEV) 6 (Kotulla et al. 2009).We adopted chemically consistent models of elliptical galaxies with a galaxy formation redshift of z = 8.We also assumed that the total initial mass of the galaxy is between 1 × 10 10 M and 5 × 10 11 M .
With these assumptions, we computed different redsequence models, which we report as horizontal lines in Fig. 7, as further outlined in the following.For all three colormagnitude plots, we show a model assuming no burst, which yields the reddest possible red sequence at the cluster redshift, as well as a model with a star formation burst at z = 2.5, exponentially declining with time.For all three clusters, the colors of the latter model are only slightly bluer than those inferred from models without a star formation burst.While these two models both nicely reproduce the colors of the reddest sources, with colors g − r 1.5, they hardly reproduce the bluer colors of BCGs.This is particularly true for the BCGs KiDS 1415 and 1433.
These considerations motivated us to perform additional models with more recent bursts of star formation, which are plotted in Fig. 7 as red horizontal lines.More recent bursts of star formation at z ∼ 1.2−1.3,z ∼ 0.45, and z ∼ 0.6−0.7 are indeed needed to better reproduce the observed g − r colors and the observed r-band magnitudes of the three KiDS 0837, 1415, and 1433 BCGs, respectively, as well.
For the BCG KiDS 0837, a burst at z = 2.5 still reproduces the BCG g − r color fairly well.However, a later burst at z ∼ 1.2−1.3appears to be favored as it implies both a bluer color and a brighter magnitude, closer to those of the BCG.Similarly to the case of KiDS 1444 discussed in C22a, we note that the adopted red-sequence models are not able to reach the low magnitude r 18.8 of the KiDS 0837 BCG.This result suggests that a faster mass growth (e.g., via mergers) than predicted by the closed-box GALEV evolutionary models would be a viable alternative for the BCG.
The color-magnitude plots and red sequence modeling together favor recent bursts of star formation for the BCGs.These results agree well with the fact that the BCGs are star forming.As further discussed in Sect.4.2 and in C22a, recently accreted pristine gas may explain the observed properties of the BCGs.For KiDS 1415, in particular, our modeling implies a remarkably recent (300 Myr) burst of star formation.This result, together with the high star formation activity and high infrared luminosity of the BCG, suggests that the source may be interacting with the nearby companion(s), which are clearly visible in Fig. 1  and CO(3→2).As pointed out in C22a for the similar case of KiDS 1220, these double-horn spectra are indicative of a low concentration of the molecular gas toward the center of the galaxy, which causes the observed flux to dim at low velocities in Fig. 3.
In the following, similarly to previous studies (Wiklind et al. 1997;Salomé et al. 2015), we attempt to link the observed double-horn spectra to the structural parameters of the gas reservoir.To do this, we followed the procedure described in Sect.5.5 of C22a and references therein, to which we refer for further details.
The observed velocity-integrated S CO(J→J−1) ∆v spectrum can be related to the H 2 gas density profile n(r) and the rotational velocity V rot (r) using Eq. ( 4) of Wiklind et al. (1997), where r is the projected radius from galaxy center, and i is the inclination with respect to the line of sight (i = 0 deg means face-on).We then parameterized the gas density profile with a Toomre (1963) disk of order n = 2, which corresponds to a flat Plummer distribution, parameterized by a characteristic disk scale d and the central gas density n 0 = 3 M H 2 2πd 2 , whose value is fixed by matching the surface integral of n(r) to the total gas mass M H 2 = 5.9 × 10 10 M .
As an alternative to the flat disk, we considered an even less centrally concentrated distribution of H 2 as described by a ring.Following Wiklind et al. (1997), we modeled a ring-like distribution as the difference of two gaseous disks with characteristic scales d 1 > d 2 .Similarly to the disk model, the characteristic density 2 ) was fixed by the surface integral.We then expressed the total rotational velocity as the quadratic sum of the contributions from the gaseous component (disk or ring), the stellar bulge, and the dark matter (DM) halo.Following Wiklind et al. (1997), we modeled the stellar bulge contribution as follows: where M = 1.86 × 10 11 M (Table 1) and r eff = 7.47 kpc are the stellar mass and half-light radius of the KiDS 1433 BCG, respectively.We estimated the latter in r band using the SourceXtractor++ (v 0.16) package7 , which is the successor to the tool SeXtractor (Bertin & Arnouts 1996).We also derived a position-dependent model of the point spread function (PSF) with the PSFEx tool8 , which we used as input in SourceXtractor++ to fit the galaxy with a Sérsic profile plus the model PSF.This analysis yielded r eff = (7.47 ± 0.40) kpc, which agrees fairly well with that inferred from the mass versus size relation of galaxies at similar redshifts (van der Wel et al. 2014), and a Sérsic index n Sérsic = 3.05.The latter is lower than the characteristic n Sérsic = 4 of elliptical galaxies, but also within the broad range of values (n Sérsic ∼ 1−10) that has been found for distant BCGs in previous studies (e.g., Bai et al. 2014;Chu et al. 2021).Following C22a, we then expressed the DM halo contribution as in Eq. ( 5) of Castignani et al. (2012), which corresponds to the universal rotation curve halo model by Salucci et al. (2007).In order to reproduce a flattening of the rotational velocity profile at a large distance r ∼ 100 kpc, as is the case for massive galaxies such as KiDS 1433 BCG, we assumed a DM mass of M DM = 4.7 × 10 12 M , which corresponds to a stellar to DM halo mass ratio of M /M DM = 4% (see, e.g., Girelli et al. 2020), and implies a flattening of the rotation velocity to a value of ∼300 km s −1 , which is typical of massive ellipticals such as KiDS 1433.In order to determine the structural parameters of the H 2 reservoir, we then varied the disk scales d, d 1 , and d 2 from 50 pc up to 20 kpc.By integrating Eq. ( 5) over all radii and over the velocity support of the CO(1→0) and CO(3→2) lines separately, we related the observed spectrum to the spectrum inferred from our modeling.The overall normalization was fixed by the height of the line.Increasing the inclination angle i has the net effect of increasing the velocity separation between the two peaks in the spectrum (Wiklind et al. 1997).As M H 2 is fixed, decreasing the disk scales d, d 1 has the net effect of increasing the characteristic density n 0 , and thus the concentration, which ultimately reduces the observed depth in the spectrum at low velocities with respect to the peaks.Given these relations, we find that a ring and a disk with an inclination of i = 45 deg both well reproduce the observed spectrum.The inferred sizes are similar for the disk and ring models, d = 7 kpc and d 1 = 5 kpc, respectively, while d 1 −d 2 = 1 kpc.With these parameters, we obtain peak H 2 -gasdensity values of ∼570 and 430 M pc −2 in the case of a disk and a ring, respectively, which are typical values found for LIRGs, such as our BCGs, along the Kennicutt-Schmidt relation.
Figure 8 displays the gas surface density profiles (a) and the rotational velocity profile (b) in the case of a d = 7 kpc disk, and similar curves are obtained in the case of a ring with the above parameters.In the bottom panels c and d we instead show the observed CO(1→0) and CO(3→2) spectra, respectively, over which we plot our model curves.The adopted set of parameters reproduces the CO(1→0) and CO(3→2) spectra well, which further supports the overall modeling.Interestingly, the star-forming KiDS 1220 BCG at z = 0.39 shows a similar double-horn profile.Recently, C22a modeled the spectrum with an extended (d = 9 kpc) and low-concentration molecular gas reservoir, similarly to the KiDS 1433 BCG.The two KiDS 1220 and KiDS 1433 BCGs also have similar stellar masses, M (2−3) × 10 11 M , gas-to-stellar mass ratios M H 2 /M (0.3−0.5), and effective radii r eff 7 kpc.Our results suggest that extended molecular gas reservoirs may be common in intermediate-z and star-forming BCGs, as indeed two (i.e., 50%) out of the four z ∼ 0.4 KiDS BCGs detected in CO by C22a and in this study exhibit double-horn CO spectra, well modeled by extended gas reservoirs.Similarly, by studying a sample of local (z < 0.017) brightest group galaxies, Olivares et al. (2022) found that 39% of their sources are dominated by extended disk and ring-like structures.These results may suggest that the star formation and gas properties of star-forming BCGs of moderately rich clusters and groups do not significantly evolve over 4 Gyr of cosmic time between z ∼ 0 and z ∼ 0.4.This agrees with recent results by Orellana-González et al. (2022), who indeed showed that the average SFR of star-forming BCGs (i.e., those with SFR > 10 M yr −1 ) does not strongly evolve with redshift, within the interval z ∼ 0−0.4.Regardless of the specific cosmic time, gas-rich BCGs are caught in a phase of significant star formation above the MS levels.As suggested by Olivares et al. (2022), during this phase, the angular momentum may prevent the cold gas from being rapidly accreted, which would explain the low concentration and the large extent of the molecular gas reservoirs.

Summary and conclusions
We have studied the molecular gas and star formation properties of three star-forming BCGs at intermediate redshifts z ∼ 0.4.We selected the sources as they exhibit significant infrared emission (Sect.2.1) from a parent sample of 684 spectroscopically confirmed z > 0.3 BCGs in the equatorial KiDS-N field (Radovich et al. 2020).Together with three similarly starforming z ∼ 0.4 KiDS BCGs drawn from the same parent sample as we studied in our recent C22a work, the targeted BCGs are among the most star-forming BCGs in KiDS.
The BCGs considered in this work have indeed line-based ([O II]) and SED-based SFR 20 M yr −1 , on average, as well as dust masses and luminosities typical of LIRGs (Sects.2.2.1 and 2.2.2).Similarly, on the basis of the position of the sources in the infrared WISE color-color plane (Jarrett et al. 2017), two BCGs are classified as starbursts and one as an intermediate disk, which confirms their significant ongoing star formation activity (Sect.2.2.3).
We investigated the stellar, star formation, and gas properties of the three BCGs with the aim of better understanding the star formation fueling and stellar mass assembly of distant BCGs.We performed an analysis of the color-magnitude plots for the cluster members (Sect.4.6).Our red-sequence modeling favors a recent burst of star formation for the BCGs.For KiDS 1415, in particular, our modeling implies a remarkably recent (300 Myr) burst of star formation.This result, together with the high star formation activity and high infrared luminosity of the BCG, suggests that the source may be interacting with its nearby companion(s).
In addition, we searched for molecular gas feeding the star formation in the three KiDS BCGs of this work.To do this, we observed them with the IRAM 30 m telescope at the first three CO transitions (Sect.3).We clearly detected double-horn CO(1→0) and CO(3→2) emission lines for one of our targets, that is, the KiDS 1433 BCG (z = 0.3546).Our observations and fits of the CO lines imply a large molecular gas reservoir of M H 2 = (5.9± 1.2) × 10 10 M and a high gas-to-stellar mass ratio M H 2 /M = (0.32 +0.12 −0.10 ).We thus increase the still limited sample of distant BCGs with detections in multiple CO transitions.The double-horn emission for KiDS 1433 implies a low gas concentration, while a modeling of the spectra yields an extended molecular gas reservoir, with a characteristic radius of ∼(5−7) kpc (Sect.4.7), which is reminiscent of the mature extended-disk phase observed in some local BCGs (Olivares et al. 2022), as well as in the star-forming KiDS 1220 BCG (z = 0.3886, C22a).These results suggest that extended molecular gas reservoirs may be common in intermediate-z and star-forming BCGs.
For the other two KiDS BCGs that we observed with the IRAM-30 m telescope in this work, we were able to set robust upper limits of M H 2 /M < 0.07 and <0.23, which are among the lowest for distant BCGs.We then combined our observations with available stellar, star formation, and dust properties of the targeted BCGs, and compared them with a sample of ∼100 distant cluster galaxies (Sect.4.1), including additional intermediate-z BCGs, with observations in CO from the literature.In particular, we confirm at higher significance than in our recent C22a work the correlation between the excitation ratio, the sSFR, and the star formation efficiency, that is, the inverse of the depletion time, separately (Sect.4.5).
In summary, we suggest that the intermediate-z KiDS BCGs that are detected in CO, which have high τ dep , M H 2 /M , and M H 2 /M dust , are observed in a peculiar and rare phase of the BCG evolution, in which a substantial amount of the H 2 gas has recently been accreted, but still not efficiently converted into stars.As has been suggested in C22a, the accretion may have occurred via interaction with the cluster companions.The other star-forming BCGs (KiDS 0837 and 1415), with upper limits in CO alone, may be at a different and possibly later stage of the BCG evolution, in which the H 2 gas reservoir has been consumed or is undergoing rapid consumption on a timescale τ dep (1−2) Gyr, while star formation is still ongoing.These results suggest that molecular gas depletion precedes star formation quenching in intermediate-z and star-forming BCGs (Sect. 4.3).This result agrees with what Castignani et al. (2022c) found locally for early-type galaxies in the field of the Virgo cluster.

Fig. 1 .
Fig. 1.Color-composite 40 × 40 KiDS images (red = i band, green = r band, and blue = g band) with OmegaCAM at the VST, centered around the BCGs KiDS 0837 (left), KiDS 1415 (center), and KiDS 1433 (right).North is up, and east is to the left.

Fig. 2 .
Fig. 2. Ultraviolet to far-infrared SEDs of the three BCGs of this work taken from GAMA DR3 (Driver et al. 2018).Dust-attenuated (red curves) and dust-unattenuated (blue curves) MAGPHYS fits are overplotted.Filled black dots and open circles are the observed and model photometry, respectively.The bottom panels show the residuals.

Fig. 3 .
Fig.3.Baseline-subtracted CO(1→0) and CO(3→2) spectra of KiDS 1433.Solid lines show the double Gaussian fits to the emission lines.See text for details.The flux (y-axis) is plotted against the relative velocity with respect to the BCG redshift, as in Table1(bottom x-axis) and the observer frame frequency (top x-axis).

Fig. 5 .
Fig. 5. SFR (top), molecular gas-to-stellar mass ratio (center), and depletion time (bottom) as a function of the stellar mass for distant BCGs and cluster galaxies observed in CO.The y-axis values are all normalized to the corresponding MS values using the relations by Speagle et al. (2014) and Tacconi et al. (2018).The horizontal dashed lines correspond to y-axis values equal to unity, and the dotted lines are the fiducial uncertainties associated with the MS.The uncertainty is chosen equal to ±0.48 dex for the top and bottom panels because the MS is commonly identified by 1/3 < SFR/SFR MS < 3.For the central panel, the plotted uncertainties are estimated at redshift z = 0.5.The color-coding for the data points is the same as in Fig. 4.
Fig. 6.Excitation ratios r J1 plotted against the sSFR (left) and the depletion time (right) for intermediate-redshift BCGs with multiple CO(J → J − 1) detections.The color-code for the BCGs and the comparison sources (ASPECS, BASS, and XCOLD GASS) is reported at the top of the two panels.

Fig. 7 .
Fig. 7. Color-magnitude plots for the sources in the field of the BCGs KiDS 0837 (top), KiDS 1415 (center), and KiDS 1433 (bottom).Field galaxies (gray points), photometrically selected cluster members (blue points), color-coded according to their cluster membership probabilities, and the BCGs (green stars) are highlighted.Red-sequence models are shown as horizontal lines.See the legend and text for further details.

Fig. 8 .
Fig. 8. Modeling of the double-horn emission lines for the KiDS 1433.(a) H 2 gas density profile.(b) Rotational velocity curve.(c,d) Observed CO(1→0) and CO(3→2) spectra and modeling.Red dots and error bars are the observed fluxes and their uncertainties, respectively.Dashed and solid black lines correspond to the disk and ring models, respectively.See the legend and text for further details.

Table 1 .
Properties of our targets.

Table 2 .
Results of our CO observations.
Table 1 (bottom x-axis) and the observer frame frequency (top x-axis).