
\section{Introduction}\label{sec:intro}
The cosmic star formation rate density increases from early times to its peak at $z\sim2$ called cosmic noon, then decreases progressively to present day \citep[for a review, ][]{2014ARA&A..52..415M}. Across cosmic history, the most massive galaxies  (stellar mass $M_\star\gtrsim10^{11}M_\odot$) have formed the bulk of their $M_\star$ around or earlier than the cosmic noon and ceased star formation in late epochs \citep[e.g.,][]{2010MNRAS.404.1775T,2013ApJ...777...18M,2014ApJ...783...85T,2015MNRAS.448.3484M,2017A&A...605A..70D}, while the formation of less massive galaxies continues for a longer period. Molecular gas is a key factor in shaping the history of galaxy assembly, as it is the immediate material of star formation \citep[for a review, ][]{2013ARA&A..51..105C}. It has been suggested that galaxies acquire molecular gas via accretion from intergalactic medium \citep[e.g.,][]{2009Natur.457..451D,2015Natur.525..496N} or mergers to fuel their star formation and central black hole growth \citep[e.g.,][]{2008ApJS..175..356H}. Observations of molecular gas content in local and distant galaxies thus provide important clues about galaxy formation and evolution.

\par
CO rotational transition lines are commonly used to trace the molecular gas content of galaxies because molecular hydrogen is a poor emitter and CO is the second most abundant molecule in the interstellar medium (ISM). Numerous observations of CO line emission have been conducted to obtain statistical samples of galaxies at intermediate and high redshifts ($z\sim1-3$) to study the relation between gas and other properties of galaxies, such as $M_\star$ and star formation rate (SFR). CO observations of high redshift galaxies typically select massive normal star-forming galaxies \citep[SFGs, e.g.,][]{2010ApJ...714L.118D,2010Natur.463..781T} or the most extreme starbursting galaxies \citep[e.g.,][]{2003ApJ...597L.113N,2010ApJ...724..233E,2013MNRAS.429.3047B}. These surveys have successfully established scaling relations that describe how galaxy properties evolve with the gas mass \citep[e.g.,][]{2015ApJ...800...20G,2018ApJ...853..179T,2020ARA&A..58..157T}. 

\par
Observations of a large sample of galaxies have revealed a tight correlation between $M_\star$ and SFR called the main sequence \citep[MS, e.g., ][]{2007A&A...468...33E,2010ApJ...714L.118D,2011ApJ...742...96W,2012ApJ...754L..29W,2014ApJS..214...15S,2015A&A...575A..74S}. SFGs on the MS dominate the cosmic star formation \citep{2011ApJ...739L..40R,2012ApJ...747L..31S}, implying the existence of mechanisms that regulate star formation from gas and the steady build-up of galaxies \citep[e.g.,][]{2013ApJ...772..119L}. Below the scatter of MS, there is another population of passive or quiescent galaxies (QGs), with little or no ongoing star formation \citep[e.g., ][]{2015ApJS..219....8C}. Many physical processes have been proposed to explain the shutdown of star formation, including starvation of gas either due to removal of cold gas by stellar and supermassive black hole feedback \citep[e.g.,][]{2006ApJS..163....1H,2010MNRAS.401....7H} or rapid gas consumption by vigorous star formation  \citep[e.g.,][]{2019MNRAS.483.5548G}, and inability of converting gas into stars due to change in the ISM condition \citep{2018ApJ...861..123F} or stabilization against collapse and fragmentation \citep[e.g.,][]{2009ApJ...707..250M}.
In order to understand the role of these processes in galaxy quenching, it is essential to characterize the gas properties such as gas fraction, depletion timescale, and kinematics.

\par
In the case of massive galaxies, this requires detecting CO emission at $z\gtrsim1$ which corresponds to the epoch when many of them are undergoing formation and subsequent quenching.  Most of the aforementioned high redshift CO observations target pre-selected MS galaxies or starbursts, thus current scaling relations only account for SFGs. Whether these relations are still valid when extrapolated to passive galaxies remains unknown. Recent targeted observations of CO emission from QGs have shown a mixed picture \citep[e.g.,][]{2021ApJ...908...54W,2021ApJ...909L..11B}. An unbiased CO survey of galaxies in various stages of star formation is needed to investigate the evolution of gas content before and after galaxy quenching. By detecting both low and medium to high-$J$ transitions, it is also possible to study the redshift evolution of molecular gas in galaxies. Several studies have used ALMA and JVLA to perform blind line searches towards deep cosmological fields or ALMA calibrators \citep{2016ApJ...833...67W,2019ApJ...882..139G,2019ApJ...872....7R,2020ApJ...896L..21R,2023MNRAS.519...34H}. However, because of the large amount of telescope time needed, the surveyed area and the number of detected sources are still very limited. 

\par
In this paper, we present serendipitous detections and physical properties of four CO emitters in the vicinity of a bright SMG SSA22-AzTEC26 from ALMA band 3 observation. In Section \ref{sec:sample} we describe the ALMA observation and data analysis. In Section \ref{sec:result} we present the extracted galaxy properties utilizing multiwavelength data in the SSA22 field. In Section \ref{sec:discussion}, we discuss the gas content of CO-selected galaxies in the context of galaxy evolution and then summarize in Section \ref{sec:summary}.

\par
Throughout this paper we adopt \citet{2003PASP..115..763C} initial mass function and Planck 2018 cosmology \citep[$H_0=67.4$ km s$^{-1}$ Mpc$^{-1}$, $\Omega_{m}=0.315$ and $\Omega_{\mathrm{\Lambda}}=0.685$, ][]{2020A&A...641A...1P}.

 \section{Sample and data analysis}\label{sec:sample}
 
\subsection{CO emitters from ALMA band 3 observation}\label{subsec:sample}

\subsubsection{ALMA data}\label{subsubsec:almaobs}
SSA22-AzTEC26 is a bright SMG first discovered by the AzTEC/ASTE 1.1 mm extragalactic survey in the SSA22 field \citep{2009Natur.459...61T,2014MNRAS.440.3462U}. The ALMA band 3 spectral scan observations (Project ID: 2019.1.01102.S, PI: Umehata)  are conducted from 20 Mar 2020 to 03 Apr 2020, towards the sky position of SSA22-AzTEC26 (RA 22:17:13.34, Dec 00:26:51.66). The observations had a maximum baseline of 313.7 m and continuously cover the sky frequency range from 84.5 GHz to 113.7 GHz with five tunings. The total integration time is 6.2 hours. J2217+0220 and J2206-0031 were observed for phase calibration. To calibrate flux and bandpass, J0006-0623, J2253+1608, J2258-2758 and J1924-2914 were observed.

\par
We use Common Astronomy Software Applications \citep[CASA, ][]{2022PASP..134k4501C} to reduce the data and perform imaging. A spectral cube is created using \texttt{tclean} task with natural weighting and a 100 km s$^{-1}$ channel width. The resulting band 3 spectral cube has $\sim1.5$  arcmin diameter field of view with a minimum beam response of 20\% of the field center. The synthesized beam size is $3.80\times2.45$ arcsec with a position angle (PA) of 57.7 degrees to $4.38\times3.65$ arcsec with PA of 81.2 degrees, depending on the sky frequencies. The RMS noise level range is $0.06-0.18$ mJy beam$^{-1}$ before primary beam correction.

\par
We also include band 4 and band 7 data (Project ID: 2021.1.01207.S, PI: Umehata) in our analysis. The band 4 observations cover the frequency ranges of 143.0-146.7 GHz and 154.9-158.7 GHz with a single tuning and total integration time of 2.5 hours. We process the data in a similar manner as band 3. The resulting band 4 cube has the synthesized beam size of $0.87\times0.79$ arcsec to $0.96\times0.88$ arcsec with PA of $\sim61$ degrees, and RMS noise level of $0.05-0.09$ mJy beam$^{-1}$ in $\sim40$ arcsec field of view. The band 7 data cover frequency range 340.4-356.2 GHz with the beam size of $0.63\times0.52$ arcsec with PA of 88.8 degrees and RMS noise of 0.09 mJy beam$^{-1}$ in the 860 $\mu$m continuum map. 


\subsubsection{Source detection and line measurements}\label{subsubsec:detection}
\par
Motivated by \citet{2016ApJ...833...67W}, we smooth the band 3 cube with a 300 km s$^{-1}$ full width at half maximum (FWHM) Gaussian filter and use the \texttt{DAOStarFinder} routine in the \texttt{photutils} package to detect all positive and negative peaks above peak signal-to-noise ratio (SNR) of 1.5 in the smoothed cube. The central AzTEC source is masked because its continuum emission is visible in the smoothed cube and contaminates source statistics. We show the number of detected peaks as a function of peak SNR in the left panel of Figure \ref{fig:fig1}. No negative peak is found above SNR$=5.2$. 
The purity of detection as a function of SNR is defined as
\begin{equation}\label{eq1}
    \mathrm{Purity} = 1-\frac{N_\mathrm{negative}}{N_\mathrm{positive}}.
\end{equation}
We fit the measured purity as a function of SNR with
\begin{equation}\label{eq2}
    P(\mathrm{SNR})=\frac{1}{2}\mathrm{erf}\left(\frac{\mathrm{SNR}-a}{b}\right)+0.5,
\end{equation} and determine the free parameters $a=4.72$ and $b=0.71$, which gives a purity of 83\% at SNR=5.2. Using this as a detection threshold, we find four sources with SNR$>5.2$ (Figure \ref{fig:fig2}). We have also tested other filter FWHM values from 100 to 500 km s$^{-1}$ but do not find any new detection.

\begin{figure*}[!htp]
    \plotone{purity.png}
    \caption{Left: The number of detections as a function of positive or negative peak SNR. The SNR$>5.2$ detection threshold is marked as vertical dotted line. Right: Purity as a function of peak SNR. The best fit relation of Equation \ref{eq2} is shown as a red curve.}
    \label{fig:fig1}
\end{figure*}

\begin{figure}[ht]
    \plotone{IRAC.png}
    \caption{IRAC 3.6 $\mu$m image of the SSA22-AzTEC26 field. The sky locations of SSA22-AzTEC26 and the four line emitters are marked with black empty circles. The large white circle indicates the ALMA band 3 field of view.}
    \label{fig:fig2}
\end{figure}

\par
For each source, we fit a single Gaussian profile to the 1D spectrum to derive the observed line frequency and FWHM (the bottom row in Figure \ref{fig:fig3}). The velocity integrated fluxes are measured by fitting a 2D Gaussian profile to the zeroth moment map, which sums along the spectral axis over a velocity range $\pm$FWHM from the line center. The measured properties are listed in Table \ref{table:table1}. Sources 1, 2, and 4 are covered by the band 4 observation, so we also report the 2 mm continuum flux densities for these three sources.  Source 1 has a second line detection at 157.548 GHz in the band 4 cube, and we will show that this is the redshifted [CI] line in the next section. 

\begin{deluxetable*}{cccccccccc}
\tablenum{1}
\caption{\label{table:table1}Result of blind line emitter search in the band 3 cube and source measurements. The 3mm fluxes and the 2 mm flux of source 1 are measured using line-free spws. For non-detection, we report $3\sigma$ upper limit. Source 1 has another line detection at 157.548 GHz, with $S\Delta v=0.320\pm0.044$ Jy km s$^{-1}$ and FWHM of $683\pm201$ km s$^{-1}$.}
\tablewidth{0pt}
    \tablehead{\colhead{ID} & \colhead{RA} & \colhead{Dec} &\colhead{$\nu_{\mathrm{obs}}$} & \colhead{SNR} & \colhead{FWHM} & \colhead{$S\Delta v$} &\colhead{$S_\mathrm{3mm}$} &\colhead{$S_\mathrm{2mm}$} & \colhead{$S_{860\mu\mathrm{m}}$}\\
    \colhead{} & \colhead{} & \colhead{} & \colhead{[GHz]} & \colhead{} & \colhead{[km s$^{-1}$]} & \colhead{[Jy km s$^{-1}$]} & \colhead{[$\mu$Jy]} & \colhead{[$\mu$Jy]} & \colhead{[mJy]}}
       
         \startdata
         1&22:17:12.36 &00:26:43.59&110.698  & 5.9&697$\pm$153 &0.569$\pm$0.075&$14.8\pm2.4$ &$90\pm18$ & -\\ 
         2&22:17:13.00 &00:26:54.06&109.151 &7.8 &635$\pm$137&0.358$\pm$0.036&$<13$ &$39\pm7.3$ & $0.815\pm0.194$\\
         3&22:17:15.01 &00:27:13.32&91.785 &7.3 &252$\pm$85&0.450$\pm$0.075& $21.7\pm3.1$ & - & -\\
         4&22:17:11.93 &00:27:14.82&107.414 &5.5 &630$\pm$123&0.713$\pm$0.141&$<27$ & $<100$ & - \\
        \enddata
\end{deluxetable*}

\subsection{Ancillary Data}\label{subsec:ancdata}
\subsubsection{Optical and NIR data}\label{subsubsec:optnir}
The SSA22 field has rich multiwavelength data coverage. We collect ground-based optical to NIR images from $u$ to $K$ bands from archives and literature. The thumbnails of the sources are shown in the top row of Figure \ref{fig:fig3}. The astrometry is corrected by stacking cutout images at the position of 1.1 mm continuum sources detected by ALMA in the SSA22 field (B. Hatsukade, private communication), and fitting a 2D circular Gaussian profile to the stacked image to determine the average positional offset between the ALMA source and its counterpart in each band. Sources 1, 2 and 4 have a coincident $K$ band counterpart with a centroid offset of 0.27, 0.23 and 0.14 arcsec, respectively, comparable with a position uncertainty of $\sim0.2$ arcsec of the CO detection. The NIR counterpart of source 3 lies 0.62 arcsec southeast of the CO position, making the association questionable.

\begin{figure*}[ht]
    \centering
    \includegraphics{stamps.png}
    \caption{
    Top row: Color-composite (red: $K$, green: $z'$, blue: $B$) cutout images of the line emitters. Each image is $10\times10$ arcsec in size. The zeroth moment maps of band 3 emission line are overlaid as white contours. The corresponding beam size is shown with white ellipses.
    Second row: the zeroth moment maps of band 3 emission line.  The zeroth moment map of band 4 157.548 GHz emission line of source 1 is shown as cyan contours. For source 2 we show the band 7 continuum image as magenta contours. All contours start from $\pm1.5\sigma$ with step of $1.5\sigma$, and the corresponding beam size is shown with the same color.
    Third row: band 3 continuum images are shown with background and purple contours. The band 4 images of sources 1 and 2 are overlaid as black contours.
    Bottom row: 1D band 3 spectra at the position of each source. The best-fit single Gaussian model is shown with black shaded area. For source 1, the band 4 emission line is shown in cyan color.  For source 4, we fit two Gaussian components and show the corresponding zeroth moment contour in the second row with the same color.}
    \label{fig:fig3}
\end{figure*}

\par
We measure the flux densities using a 2 arcsec diameter circular aperture placed at the CO position in bands shorter than IRAC. As seen in Figure \ref{fig:fig3}, the flux of source 4 is contaminated by two blue neighboring sources. Therefore, we model the image using three Sersic profiles simultaneously and subtract the model of neighbors before flux measurement as shown in Figure \ref{fig:fig4}. Then we aperture-correct the fluxes using the point spread function (PSF) and correct for galactic extinction using \citet{1998ApJ...500..525S} dust map and \citet{1989ApJ...345..245C} Milky way extinction curve.

\begin{figure}
    \plotone{gal4deblend.png}
    \caption{Examples of neighbor subtraction of galaxy 4. The left, middle, and right columns show original data, best-fit Sersic model, and neighbor-subtracted data, respectively. The 2 arcsec diameter aperture used to measure the flux is shown as red circle.}
    \label{fig:fig4}
\end{figure}

\subsubsection{Infrared data}\label{subsubsec:irdeblending}
 For long wavelength NIR (3.6   $\mu$m) to mid-infrared (MIR), we obtain IRAC 3.6-8 $\mu$m and MIPS 24$\mu$m images from the Spizer archive. At far-infrared (FIR) wavelengths, we use Herschel SPIRE 250-500 $\mu$m images from the ESA Herschel archive. Flux densities are derived from PSF photometry to overcome source blending issues. For the IRAC data, we apply the method introduced by \citet{2012ApJS..203...23H} which iteratively decomposites the image into the combination of ideal point sources and background noise. We consider extracted fluxes within 1 arcsec (1.5 pixels of IRAC images) as associated with the CO emitter. For Spitzer/MIPS and Herschel/SPIRE images, we measure the flux densities by directly fitting PSF models to the image with IRAC 4.5 $\mu$m and 1.1 mm source positions as prior, following the method of \citet{2017MNRAS.464..885H}. Since our sources lie on the edge of the MIPS coverage, we also consider sources detected at SNR$>5$ in IRAC 8 $\mu$m image when deblending the SPIRE images to add dusty sources that are not covered by MIPS but might appear in SPIRE images.  The measured fluxes are listed in Table \ref{table:table2}.

\begin{deluxetable*}{ccccccc}[ht]
 \tablenum{2}
\caption{\label{table:table2}Multiwavelength photometry of the four CO emitters in unit of $\mu$Jy, extracted at the CO position.}
\tablewidth{0pt}
\tablehead{\colhead{Band} & \colhead{Instrument} &\colhead{1} &\colhead{2} &\colhead{3} &\colhead{4} &\colhead{Reference}}
\startdata
$u$&CFHT/MegaCam &-&$0.125\pm0.020$&$1.183\pm0.019$&$0.076\pm0.019$&1 \\
$B$&Subaru/SuprimeCam &$0.057\pm0.006$&$0.113\pm0.007$&$1.388\pm0.007$&$0.019\pm0.007$&2 \\
$V$&... &$0.109\pm0.009$&$0.199\pm0.010$&$2.303\pm0.009$&$0.091\pm0.009$&2 \\
$R$&... &$0.142\pm0.014$&$0.371\pm0.013$&$3.316\pm0.013$&$0.333\pm0.015$&3 \\
$i'$&... &$0.186\pm0.018$&$0.748\pm0.018$&$5.853\pm0.020$&$0.885\pm0.021$&2 \\
$z'$&... &$0.311\pm0.034$&$1.609\pm0.036$&$8.122\pm0.038$&$2.679\pm0.037$&2 \\
NB359&... &$0.067\pm0.019$&$0.058\pm0.020$&$1.223\pm0.019$&$0.116\pm0.019$&4 \\
NB497&... &$0.099\pm0.014$&$0.124\pm0.015$&$1.874\pm0.015$&$0.068\pm0.016$&5 \\
NB816&... &$0.249\pm0.036$&$1.023\pm0.037$&$6.455\pm0.036$&$1.322\pm0.039$&2 \\
NB912&... &$0.265\pm0.029$&$1.462\pm0.030$&$7.906\pm0.029$&$2.977\pm0.034$&2 \\
$g$&Subaru/HSC &$0.083\pm0.017$&$0.133\pm0.016$&$1.828\pm0.017$&$0.088\pm0.017$&6 \\
$r$&... &$0.115\pm0.031$&$0.280\pm0.033$&$3.077\pm0.034$&$0.303\pm0.035$&6 \\
$i$&... &$0.177\pm0.027$&$0.793\pm0.030$&$6.245\pm0.030$&$0.948\pm0.029$&6 \\
$z$&... &$0.284\pm0.079$&$1.668\pm0.076$&$7.797\pm0.077$&$2.677\pm0.079$&6 \\
$Y$&... &-&$2.123\pm0.126$&$8.042\pm0.123$&$4.182\pm0.122$&6 \\
$J$&UKIRT/WFCAM &$1.019\pm0.206$&$4.137\pm0.214$&$13.052\pm0.214$&$8.515\pm0.210$&7 \\
$K$&... &$4.071\pm0.322$&$13.838\pm0.310$&$25.458\pm0.340$&$32.593\pm0.328$&7 \\
IRAC1&Spitzer/IRAC &$11.565\pm0.220$&$25.925\pm0.192$&$44.430\pm0.704$&$67.229\pm0.556$&8,9 \\
IRAC2&... &$18.926\pm0.851$&$25.833\pm1.121$&$44.152\pm0.815$&$59.478\pm2.734$&8,9\\
IRAC3&... &$17.952\pm2.103$&$17.459\pm1.747$&$40.255\pm1.980$&$37.444\pm2.027$&8,9 \\
IRAC4&... &$15.686\pm2.269$&$18.327\pm2.857$&$28.632\pm2.483$&$32.981\pm3.147$&8,9\\
MIPS1&Spitzer/MIPS &$305\pm149$&$20\pm44$&$392\pm169$&$115\pm189$&8,9 \\
PSW&Herschel/SPIRE &$6133\pm4584$&$4931\pm4861$&$4997\pm5191$&$3432\pm4442$&10,11 \\
PMW&... &$6437\pm7227$&$11129\pm6856$&$5741\pm6574$&$2406\pm4033$&10,11 \\
PLW&... &$1941\pm3710$&$4385\pm5905$&$2523\pm4818$&$1284\pm2745$&10,11 \\
\enddata
\tablerefs{
1: K. Mawatari, private communication;
2: \citet{2011MNRAS.412.2579N};
3: \citet{2004AJ....128.2073H};
4: \citet{2009ApJ...692.1287I};
5: \citet{2012AJ....143...79Y};
6: \citet{2018PASJ...70S...8A} 
7: \citet{2007MNRAS.379.1599L}; 
8: \citet{2009ApJ...692.1561W};
9: \url{https://sha.ipac.caltech.edu/applications/Spitzer/SHA/};
10: \citet{2016MNRAS.460.3861K};
11: \url{http://herschel.esac.esa.int/Science_Archive.shtml}
}
\tablecomments{``..." means the same as the previous band.}
\end{deluxetable*}

\subsection{SED modelling}\label{subsec:SEDfitting}
We adopt the methodology of \texttt{prospector} \citep{2019ApJ...877..140L} to model the observed SED, which fits all physical parameters simultaneously and yields the joint probability density function (PDF) as result. We model the dust attenuation for old and young stars with two \citet{2009A&A...507.1793N} dust attenuation curves, which are \citet{2000ApJ...533..682C} dust extinction curve multiplied by a powerlaw modification $(\lambda_\mathrm{rest}/550\mathrm{nm})^{\beta}$. For the two stellar populations separated by age, the slope $\beta$ and amplitude of each curve can vary independently within a flat prior [-0.7, 0.4] for $\beta$ and [0, 10] for $E(B-V)$, while the dust attenuation of old stars is forced to be smaller than that of young stars. This adds some additional flexibility to take possible complex dust geometry into account. To implement the modification, SED models are built with the \texttt{CIGALE} code \citep{2019A&A...622A.103B}. Before line identification, we let redshift vary freely within a flat prior [0, 8] (i.e., photometric redshift mode of SED fitting).

\section{Results}\label{sec:result}
\subsection{Line identification}\label{subsec:zsolve}
\begin{figure*}[htp!]
    \plotone{zsolve.png}
    \caption{Observed frequency as a function of redshift for different spectral lines. CO lines are shown with colored curves and [CI] line is shown with dashed black curves. Blue shaded areas show the frequency coverage of the ALMA band 3 observation. We mark the observed frequency of each galaxy with a horizontal black line and the corresponding redshift of each solution with a vertical dotted line. The PDFs of photometric redshift are illustrated by lines and shaded areas, with black color for the SED fitting method described in Section \ref{subsec:SEDfitting} and orange color for \texttt{EASY}. A successful solution requires the vertical line to locate within the photometric redshift PDF.}
    \label{fig:fig5}
\end{figure*}

Each of all sources has only one robust line detection in the band 3 cube. The most probable identification of the detected line will be one of CO $J_\mathrm{up}=4,3,2,1$ from galaxies at $z\lesssim3-4$ because our $\sim29$ GHz bandwidth will cover two CO lines or one CO and [CI]($^3\mathrm{P}_1-^3\mathrm{P}_0$) atomic carbon emission line for $z\gtrsim3-4$. To find the redshift solution, we plot the marginalized PDF of photometric redshift from our SED fitting and \texttt{EASY}\footnote{\url{https://github.com/gbrammer/eazy-py}} code \citep{2008ApJ...686.1503B}, together with observed frequency as a function of redshift for CO lines and [CI]($^3\mathrm{P}_1-^3\mathrm{P}_0$) line in Figure \ref{fig:fig5}.

\par
Based on the result of SED fitting, we find that all four sources are distant galaxies at a cosmological distance.  For galaxies 1, 2 and 4, there is a unique solution of CO line redshift that is allowed by the photometric redshift PDF. The redshift of 2.124 given by CO J=3-2 line of galaxy 1 agrees with the second line in band 4 being redshifted [CI]($^3\mathrm{P}_1-^3\mathrm{P}_0$), which yields consistent redshift $z_{\mathrm{[CI]}}=2.124$.  The photometric redshift PDF of galaxy 3 disfavors all solutions. Since the flux is extracted at the CO position, we check the result by changing the position of the photometry aperture to the center of the nearby NIR object, then the photometric redshift of $0.63^{+0.07}_{-0.06}$ still does not overlap with any solution. Thus, the flux is likely to be dominated by a NIR source that is not associated with the CO emission. Given only one line detection within the $\sim29$ GHz bandwidth, a possible solution will be 
$z_{\mathrm{CO(2-1)}}=1.512$ or $z_\mathrm{CO(3-2)}=2.767$ but need to be verified by detecting other spectral lines in future observation \citep[e.g.,][]{2014ApJ...781L..39T,2021ApJ...917...94M}.

\begin{deluxetable}{ccccc}
\tablenum{3}
\caption{\label{table:table3}Results of line identification.}
\tablehead{\colhead{ID} &\colhead{CO transition}&\colhead{$z_\mathrm{CO}$}  &\colhead{$L'_{\mathrm{CO}(J_\mathrm{up}-J_\mathrm{low})}$}  &\colhead{$M_\mathrm{mol}$\tablenotemark{1}} \\
\colhead{} &\colhead{[$J_\mathrm{up}-J_\mathrm{low}$]} &\colhead{} &\colhead{[$10^{10}$ K km s$^{-1}$ pc$^2$]} &\colhead{[$10^{10}M_\odot$]}}
\startdata
        1 &3-2 &2.124 &$1.45\pm0.19$&$12.5\pm1.6$\\
        2 &2-1 &1.113 &$0.62\pm0.06$&$2.9\pm0.3$\\
       4 &2-1 &1.146 &$1.31\pm0.26$&$6.2\pm1.2$\\
\enddata
\tablenotetext{1}{CO-based molecular gas mass (Section \ref{subsec:mgas}). Alternatively, galaxy 1 has  $M_\mathrm{mol}=(9.4\pm1.3)\times10^{10}M_\odot$ based on its [CI](1-0) line luminosity $L'_\mathrm{[CI](1-0)}=(4.04\pm0.55)\times10^9$ K km s$^{-1}$ pc$^2$.}
\end{deluxetable}

\par
The line luminosity is calculated from the velocity integrated line flux $S\Delta v$, following \citet{2005ARA&A..43..677S}:
\begin{equation}\label{eq3}
\frac{L'}{\text{K km s}^{-1}\text{ pc}^2} = 
\frac{3.25\times10^7S\Delta v}{\text{Jy km s}^{-1}}
\frac{D_\mathrm{L}^2}{(1+z)^3\nu_\mathrm{obs}^2}\left[\frac{\mathrm{GHz}}{\mathrm{Mpc}}\right]^2,
\end{equation}
where $D_\mathrm{L}$ is the luminosity distance in Mpc corresponding to redshift $z$, and $\nu_\mathrm{obs}$ is the observed frequency of the line in GHz. We list the result of line identification in Table \ref{table:table3}.  With the spectroscopic redshift determined from the CO line, we fix the redshift and refit the photometry for galaxies 1, 2 and 4, including band 3, 4 and 7 continuum measurements. In the following sections, we will focus on the derived physical properties of these three galaxies. 

\subsection{Stellar mass and SFR from SED modelling}\label{subsec:mstarsfr}
\begin{deluxetable}{cccccc}
\tablenum{4}
\caption{\label{table:table4}Derived physical properties from SED fitting. The SFRs are averaged over the past 30 Myrs.} 
\tablehead{\colhead{ID} &\colhead{$M_\star$} &\colhead{SFR} &\colhead{$M_\mathrm{dust}$}&\colhead{$L_\mathrm{dust}$} &\colhead{$A_V$}\\
\colhead{}&\colhead{[$10^{10}M_\odot$]}&\colhead{[$M_\odot$ yr$^{-1}$]}&\colhead{[$10^{8}M_\odot$]}&\colhead{[$10^{11}L_\odot$]}&\colhead{[mag]}}
\startdata
        1 &$9.91^{+2.22}_{-2.54}$ &$100^{+48}_{-36}$ &$3.21^{+0.64}_{-0.56}$ &$12.11^{+3.57}_{-3.55}$&$2.49^{+0.24}_{-0.31}$\\
        2 &$4.89^{+0.74}_{-0.82}$&$24^{+8}_{-6}$&$2.30^{+0.79}_{-0.50}$&$3.21^{+0.80}_{-0.50}$&$2.15^{+0.19}_{-0.20}$\\
        4 &$22.76^{+2.66}_{-2.57}$&$0.2^{+0.7}_{-0.2}$ &$0.32^{+0.23}_{-0.09}$&$1.29^{+0.18}_{-0.20}$&$1.30^{+0.11}_{-0.14}$\\
\enddata
\end{deluxetable}

\begin{figure}
    \plotone{sfms.png}
    \caption{The position of the three galaxies with CO redshift in SFR-$M_\star$ plane. The $M_\star$ and SFRs are derived from SED modeling (Table. \ref{table:table4}). Circles mark galaxies detected with  CO(2-1) while squares mark CO(3-2). MS relation at $z=1.1$ and $z=2.1$ are shown with dotted and dashed lines, respectively. The shaded area indicates the scatter of MS.}
    \label{fig:fig6}
\end{figure}

\par
We show the best-fit SED and joint posterior distributions of derived parameters in Figure \ref{fig:figa1} and \ref{fig:figa2}. The results from SED fitting are summarized in Table \ref{table:table4}. All reported values and errors represent median values and 68 percent confidence intervals from marginalized posterior PDF, respectively.

\par
In Figure \ref{fig:fig6} we plot the three galaxies in the SFR-$M_\star$ plane, together with the MS relation at $z=1.1$ and $z=2.1$ from \citet{2014ApJS..214...15S}. For comparison, we also show the CO-selected galaxy sample from ASPECS-LP \citep{2019ApJ...882..140B} that surveyed a larger area at a similar depth. The ASPECS-LP $M_\star$ and SFR are derived using \texttt{MAGPHYS} code \citep{2015ApJ...806..110D}. Different SFH and dust attenuation models can lead to a $\sim0.3$ dex systematic offset in the SED-derived $M_\star$ and SFR \citep{2019A&A...632A..79B}. Tests on simulated galaxies show the non-parametric SFH used in this study gives a larger $M_\star$ than simple exponential SFH models \citep{2019ApJ...876....3L,2020ApJ...904...33L}. However, the results should still agree within error for the same object as long as each method allows a wide range of SFHs and dust curves.

\par
All three galaxies with CO-based redshift are massive, with $M_\star>10^{10.5}M_\odot$ but different levels of star formation. Galaxies 1 and 2 lie within the scatter of the MS at their redshifts. In contrast, galaxy 4 is quiescent with specific SFR $\mathrm{sSFR}=\mathrm{SFR}/M_\star=0.9_{-0.8}^{+2.8}\times10^{-12}$ yr$^{-1}$ and log($\Delta$MS)=log(SFR/SFR$_\mathrm{MS}$) = $-2.7^{+0.6}_{-0.9}$. Compared with the ASPECS-LP galaxies, galaxies 1 and 2 have comparable $M_\star$ and SFR, but galaxy 4 is a significant outlier, with SFR much lower than the median value (30 $M_\odot$ yr$^{-1}$) of galaxies detected with CO(2-1) in the ASPECS-LP survey.

\begin{figure}
    \plotone{sfr.png}
    \caption{SFR$_\mathrm{SED}$ versus SFR$_\mathrm{UV+IR}$. Dashed line denotes 1:1 equality.}
    \label{fig:fig7}
\end{figure}

\par 
We further compare our SED-based SFR with those derived from UV+IR conversions in Figure \ref{fig:fig7}. The UV+IR SFR is calculated following \citet{2005ApJ...625...23B}:
\begin{equation}\label{eq4}
    \frac{\mathrm{SFR}_\mathrm{UV+IR}}{M_\odot \text{ yr}^{-1}}=\frac{1.09\times10^{-10}(2.2L_\mathrm{UV}+L_{\mathrm{IR}}
    )}{L_\odot},
\end{equation}
where $L_\mathrm{UV}$ is the rest-frame luminosity integrated in the rest-frame wavelength range 121.6-300 nm without dust correction and $L_\mathrm{IR}$ is the luminosity integrated in the rest-frame wavelength range 8-1000 $\mu$m. Galaxies 1 and 2 have SFR$_\mathrm{SED}/$SFR$_\mathrm{UV+IR}$ of 79\% and 69\%, respectively. Comparatively, galaxy 4's SFR$_\mathrm{SED}$ is only $\sim1$\% of SFR$_\mathrm{UV+IR}$. This trend decreases while remaining visible if we change the average interval of SFR$_\mathrm{SED}$ from 30 Myrs to 100 Myrs. The difference between the two SFR methods increases with decreasing sSFR because of the contribution of dust heating from old stars \citep{2014MNRAS.445.1598H,2019ApJ...877..140L}, and this behavior is also found when using \texttt{MAGPHYS} code \citep{2019ApJ...882...65M}. Since dust attenuation of old and young stars are modeled separately and both low and high attenuation values are allowed for the two stellar populations, the error range of SFR$_\mathrm{SED}$ serves as a conservative estimate of SFR levels.  Therefore, we continue to choose SFR$_\mathrm{SED}$ as our SFR indicator. 


\subsection{Molecular gas mass}\label{subsec:mgas}
In order to derive their molecular gas mass, we first convert the CO luminosity to the ground transition assuming line ratios $r_{21}=0.76\pm0.09$ and $r_{31}=0.42\pm0.07$ appropriate for high redshift normal SFGs \citep{2015A&A...577A..46D,2016ApJ...833...70D}. The total molecular gas mass $M_\mathrm{mol}$ is calculated using the conversion factor of normal star-forming galaxies: $\alpha_{\mathrm{CO}}=M_\mathrm{mol}/L'_{\mathrm{CO(1-0)}}=$ 3.6 $M_\odot/(\mathrm{K}$ km s$^{-1}$ pc$^2)$ \citep{2010ApJ...713..686D}. The resulting gas masses are listed in Table. \ref{table:table2}. For galaxy 1 we can also estimate $M_\mathrm{mol}$ from the [CI] line flux. The conversion factor $\alpha_{\mathrm{[CI]}}=M_\mathrm{mol}/L'_{\mathrm{[CI]}}$ is calculated using the equation 11 of \citet{2022MNRAS.517..962D}:
\begin{equation}\label{eq5}
 \frac{\alpha_{\mathrm{[CI]}}}{M_\odot(\text{K km s}^{-1} \text{pc}^2)^{-1}}
 =16.8\left(\frac{X_{\mathrm{CI}}}{1.6\times10^{-5}}\right)^{-1}\left(\frac{Q_{10}}{0.48}\right)^{-1},
\end{equation}
where $X_\mathrm{CI}$ is the average abundance ratio of atomic carbon and $Q_{10}$ is the excitation term. We adopt $X_\mathrm{CI}=10^{-4.8}$ and $Q_{10}=0.35$ following \citet{2021ApJ...909..181L}.
Using the band 4 data we find the velocity integrated [CI](1-0) line flux $S_\mathrm{[CI](1-0)}\Delta v=0.320\pm0.044$ Jy km s$^{-1}$. This translates to a [CI](1-0)-based molecular gas mass of $(9.4\pm1.3)\times10^{10}M_\odot$, roughly consistent with the CO-based gas mass.

\par
We note that molecular gas mass estimates suffer from substantial uncertainties in $\alpha_\mathrm{CO}$ and excitation.  In this study, we have adopted $\alpha_\mathrm{CO}$ and line ratios of normal SFGs. Some recent ALMA studies have revealed compact star formation with intense starburst-like conditions in high redshift MS galaxies, or ``hidden starburst in the MS'' \citep[e.g., ][]{2021MNRAS.508.5217P}. In this case, starburst-like conversions should be used to derive gas mass despite their SFR being within the scatter of the MS. If we instead assume SMG-like $\alpha_\mathrm{CO}=1.36$ $M_\odot/({\mathrm{K}}$ km s$^{-1}$ pc$^2)$, $r_{21}=0.84\pm0.13$ and $r_{31}=0.52\pm0.09$ \citep{2013MNRAS.429.3047B}, the gas mass will decrease by 2.96 and 3.28 times for CO(2-1) and CO(3-2), respectively. For [CI](1-0) line, with SMG-like $X_\mathrm{[CI]}=10^{-4.2}$ and $Q_{10}=0.48$ \citep{2011ApJ...730...18W,2018ApJ...869...27V} the estimated molecular gas mass of galaxy 1 becomes $(1.7\pm0.2)\times10^{10}M_\odot$, or 45\% of the CO-based gas mass.

\par
The uncertainties in gas mass estimates are further complicated by the discrepancy between different tracers, especially CO and dust. In the case of galaxy 1, normal SFG conversions give a gas-to-dust mass ratio $\delta_\mathrm{GDR}=387^{+102}_{-80}$, more than three times higher than $\delta_\mathrm{GDR}=110$ expected by the  $M_\star$-metallicity-redshift and metallicity-$\delta_\mathrm{GDR}$ scaling relations \citep[e.g., ][]{2008A&A...488..463M,2011ApJ...737...12L,2015ApJ...800...20G}, while SMG-like conversions bring $\delta_\mathrm{GDR}=118^{+30}_{-25}$ closer to the scaling relations. On the other hand, without changing to SMG-like factors, galaxy 2 already has $\delta_\mathrm{GDR}=126^{+40}_{-32}$ that is consistent with $\delta_\mathrm{GDR}=102$ from scaling relations. Thus, the ``hidden starburst in the MS'' scenario cannot be applied to all galaxies. Moreover, the position of galaxy 1 with respect to the $L_\mathrm{dust}-L^{'}_{\mathrm{CO(1-0)}}$ relation does not support a starburst picture. \citet{2010MNRAS.407.2091G} find $\log(L_\mathrm{dust})=1.15\log(L^{'}_{\mathrm{CO(1-0)}})+0.02\pm1.1$ for normal SFG and the same slope but intercept of $0.63\pm0.12$ for starbursts. For galaxy 1 we derive $\log(L_\mathrm{dust}/L^{'}_{\mathrm{CO(1-0)}}{}^{1.15})=-0.03^{+0.14}_{-0.16}$ under normal SFG conversion or $0.07^{+0.14}_{-0.16}$ under SMG conversion, thus it is likely that the former is more appropriate. We consider the dust-based $M_\mathrm{mol}$ as the least robust here because of low SNR ($<2$) in all Herschel bands and the uncertainty in dust SED fitting \citep[e.g.,][]{2016A&A...587A..73B}. In the following discussion, we keep using normal SFG conversions.


\section{Discussion}\label{sec:discussion}

\begin{figure*}
    \plotone{z_gas.png}
    \caption{Left: Molecular gas fraction as a function of redshift. For comparison we show the scaling relation and its scatter from \citet{2018ApJ...853..179T}, local star-forming galaxies \citep[xCOLD GASS, ][]{2017ApJS..233...22S}, high redshift star-forming galaxies \citep[PHIBSS1/2, ][]{2013ApJ...768...74T,2019A&A...622A.105F}, local QGs \citep[MASSIVE, ][]{2019MNRAS.486.1404D}, intermediate and high redshift QGs \citep{2018ApJ...860..103S,2021ApJ...908...54W,2021ApJ...909L..11B}. Triangles indicate upper limits. Right: SFR as a function of $M_\mathrm{mol}$. Symbols as the left panel.}
    \label{fig:fig8}
\end{figure*}

\subsection{Comparison with scaling relation}\label{subsec:gasscaling}
The left panel of Figure \ref{fig:fig8} shows the gas ratio $\mu_\mathrm{mol}=M_\mathrm{mol}/M_\star$ of our sample and values of SFGs and QGs across cosmic time, color-coded by $\log(M_\star/M_\odot)$. The SFG sample is collected from \citet{2018ApJ...853..179T} which uses \citet{2015ApJ...800...20G} metallicity-dependent prescription ($\alpha_\mathrm{CO}=3.8-12.4$ $M_\odot$(K km s$^{-1}$ pc$^{2}$)$^{-1}$ for $z=0-4$ and $\log(M_\star/M_\odot)=9-12$). The $M_\mathrm{mol}$ of QGs in the comparison sample is derived using Galactic $\alpha_\mathrm{CO}$. At $z=1-1.5$, the SFR$>30$ $M_\odot$ yr$^{-1}$ and $M_\star>2.5\times10^{10}M_\odot$ cuts of PHIBSS1 survey makes their pre-selection slightly higher than MS SFR at this redshift,  resulting in median $\mu_\mathrm{mol}$ of 0.79. This is higher than $\mu_\mathrm{mol}=0.61_{-0.10}^{+0.13}$ of galaxy 2 or $0.27_{-0.06}^{+0.07}$ of galaxy 4. Similar differences have been found by ASPECS-Pilot survey \citep{2016ApJ...833...70D} which found that CO emitters from their ALMA band 3 line search in the HUDF field have $\sim2\times$ lower median gas ratio than PHIBSS1. Also, in subsequent larger ASPECS-LP sample \citep{2019ApJ...882..136A}, the galaxies detected with CO(2-1) have median $\mu_\mathrm{mol}=0.48$. At $z>2$, the median gas ratio of PHIBSS1 galaxies is 1.26, comparable with $\mu_\mathrm{mol}=1.27_{-0.28}^{+0.47}$ of galaxy 1, as at this redshift their SFR and $M_\star$ cuts do not bias the sample to galaxies above the MS. Our sample further contains one galaxy that is significantly below the MS, showing a blind line search is able to cover a wider range of star formation levels.

\par
The right panel of Figure \ref{fig:fig8} shows SFR as a function of $M_\mathrm{mol}$. The two quantities are closely correlated with modest time evolution of depletion time $t_\mathrm{dep}=M_\mathrm{mol}/\mathrm{SFR}\propto(1+z)^{-0.63}$ \citep{2018ApJ...853..179T}. In our sample, galaxies 1 and 2 show a lower SFR at given $M_\mathrm{mol}$ compared to targeted SFGs in the PHIBSS1 survey, but their $t_\mathrm{dep}=1.2_{-0.4}^{+0.7}$ and $1.2_{-0.3}^{+0.4}$ Gyr agree with the scaling relation, which expects $\sim0.8$ and $\sim0.9$ Gyr, respectively. This might be caused by the difference in methods for deriving SFR and $M_\star$. As mentioned in Section \ref{subsec:mstarsfr}, the full UV-to-FIR SED fitting method used in this study typically gives lower SFR and larger $M_\star$. Such a trend is also found by ASPECS-LP \citep{2019ApJ...882..136A}, as they use \texttt{MAGPHYS} that also includes dedicated treatment of SFH and dust attenuation. In summary, while differs in target selection and methodology of analysis, our sample shows that the current scaling relations successfully describe the evolution of gas content  as a function of redshift, mass, and SFR in star-forming galaxies.

\subsection{Large gas reservoir in a massive galaxy below the MS}\label{subsec:gal4}

The most surprising finding from the analysis above is the large gas reservoir in massive galaxy 4 with low sSFR. CO/dust observations of QGs outside the local universe based on pre-selection of sSFR typically find low gas fraction of $\mu_\mathrm{gas}\lesssim10\%$ \citep[][ see also the left panel of Figure \ref{fig:fig8} for CO samples]{2015ApJ...806L..20S,2018ApJ...860..103S,2019ApJ...873L..19B,2021ApJ...908...54W,2021ApJ...910L...7C}. On the other hand, higher values of $\mu_\mathrm{mol}$ in QGs have also been reported \citep{2017ApJ...849...27R,2018ApJ...856..118H,2018NatAs...2..239G,2021ApJ...909L..11B}. Here galaxy 4 has $M_\mathrm{mol}=(6.2\pm1.2)\times10^{10}$ $M_\odot$ and $\mu_\mathrm{gas}$ of $0.27_{-0.06}^{+0.07}(\alpha_\mathrm{CO}/3.6$ $M_\odot$ $({\mathrm{K}}$ km s$^{-1}$ pc$^2)^{-1})((0.76\pm0.09)/r_{21})$, adding one more candidate to gas-rich QGs. With the results, we infer a very low star formation efficiency (SFE) of galaxy 4 with SFE$=1/t_\mathrm{dep}=0.34^{+1.10}_{-0.31}\times10^{-11}$ yr$^{-1}$. The $t_\mathrm{dep}$ is not only longer than targeted SFGs in PHIBSS surveys  but also any blindly detected galaxies from ASPECS-LP, which has the lowest SFE$=2.2^{+0.5}_{-0.3}\times10^{-10}$ yr$^{-1}$. In contrast, pre-selected QGs detected in CO at $z\gtrsim1$ show SFE closer to normal SFGs (the right panel of Figure \ref{fig:fig8})

\par
Why is galaxy 4 so inefficient in forming new stars while hosting a large gas reservoir? If the gas reservoir is leftover after the end of the main star formation episode without being depleted or destroyed, there must be some mechanisms to halt new star formation. In bulge-dominated systems, the gas disk can be stabilized against collapse \citep{2009ApJ...707..250M,2014ApJ...785...75G} so further star formation is dynamically suppressed. However, the gas fraction of galaxy 4 is too high for morphological quenching to effectively reduce star formation \citep[$\lesssim5\%$, ][]{2013MNRAS.432.1914M,2020MNRAS.495..199G}. From the SFH (Figure \ref{fig:figa2}) we see the bulk of stars were formed $\gtrsim2$ Gyr ago and the SFR drops at $z\sim2$, with mass-weighted stellar age of $3.2^{+0.5}_{-0.3}$ Gyr. Other sources such as supernova or AGN \citep{2010A&A...521A..65N} can also inject kinematic energy to the ISM and prevent collapse, but these processes take place on shorter timescales \citep[$\lesssim100$ Myr, ][]{2015A&A...574A..32G} and may also remove the cold gas, so they are at least not the main reason of galaxy 4's quiescence while keeping the high gas fraction. Finally, the $\gtrsim1$ Gyr lookback time is much longer than the remaining lifetime of the gas reservoir after quenching inferred from post starburst galaxies \citep[$\sim150$ Myr, ][]{2022ApJ...925..153B}. Thus, the gas reservoir in galaxy 4 is unlikely the remnant of the past major star formation event, but the listed mechanisms might have played a role in keeping the low SFE.

\par
Alternatively, the gas reservoir might be acquired in late times \citep{2022ApJ...940...39W}. Since galaxy 4 is already very massive ($M_\star>10^{11}M_\odot$), its halo mass will be $M_\mathrm{halo}>10^{12.5}M_\odot$ assuming a stellar-to-halo mass ratio of -1.5 dex \citep{2013ApJ...770...57B}. At this $M_\mathrm{halo}$, shock heating can suppress further cold gas supply via accretion from the surrounding environment \citep{2005MNRAS.363....2K,2006MNRAS.368....2D}. Rather than accretion, several minor mergers could have added gas to galaxy 4. According to the MS gas scaling relation, a normal galaxy with $1/100-1/4$ of galaxy 4's $M_\star$ and redshift $1.146-2$ carries $M_\mathrm{gas}\sim0.3-4.8\times10^{10}$ $M_\odot$. Integrating the merger rate found by Illustris simulation \citep{2015MNRAS.449...49R} over this mass and redshift range, there are 2.9 mergers expected and $3.4\times10^{10}$ $M_\odot$ of added molecular gas, which is only half of the gas mass of galaxy 4. Minor merger scenario also struggles to explain the low SFR, though they do not necessarily elevate SFR immediately due to dynamical effects \citep{2015MNRAS.449.3503D,2018MNRAS.476..122V}. If minor mergers really happened to galaxy 4, there should be some imprints left in it, such as stellar-gas misalignment \citep[e.g.,][]{2021ApJS..254...27K}. The CO data show a rotating gas disk might exist in galaxy 4 (Figure \ref{fig:fig3}), but currently we cannot perform kinematic diagnostics due to a lack of optical/near-infrared spectroscopic data. We conclude that none of the mechanisms discussed above can explain galaxy 4's high gas fraction alone, and future high-resolution observations are needed to reveal the nature of galaxy 4. Such observations may include high-resolution optical/NIR imaging to study stellar morphology; confirmations of quiescence and relation with neighbor sources using optical/NIR spectroscopy;  and interferometric observation of CO line kinematics, especially for diagnostics of stability, temperature, and density of its gas disk.  

\subsection{Caveats of galaxy 4's SFR estimate}\label{subsec:caveat}
\par
We caution that the current identification of galaxy 4's as a QG is only tentative.  The assertion of the low SFR is based on very faint rest-frame UV fluxes, which depend on deblending optical images to remove contamination from the two blue neighboring sources (Figure \ref{fig:fig4}). Without deblending, the SFR$_\mathrm{SED}$ will be $9_{-8}^{+15}$ $M_\odot$ yr$^{-1}$. While their photometric redshifts of $0.89_{-0.14}^{+0.07}$ and $0.29^{+0.17}_{-0.26}$ suggest the two neighbors are irrelevant to galaxy 4, it is still possible that they are actually physically associated with galaxy 4 because of uncertainties in deblended fluxes. If they are indeed at the same redshift, they would be blue dwarf galaxies with $\lesssim1\%$ of $M_\star$ of galaxy 4 and SFR of $\lesssim20$ $M_\odot$ yr$^{-1}$ in total.

\par
Besides the uncertainty in unobscured SFR traced by rest-frame UV fluxes, the obscured SFR of galaxy 4 is left unconstrained. In this $M_\star$ range, almost all of the star formation is obscured \citep[$>90\%$, e.g., ][]{2017ApJ...850..208W}, so FIR/radio observations are needed to measure the total SFR. The galaxy is not detected in Herschel SPIRE bands, which have a$1\sigma$ depth of $\sim7$ mJy. The $3\sigma$ upper limits of ALMA band 3/4 continuum fluxes are still 1-2 orders of magnitudes higher than that predicted by SED fitting (Figure \ref{fig:figa1}). The deepest radio data is from JVLA 3 GHz observations \citep{2017ApJ...850..178A} of the SSA22 field with a primary beam response of $\sim0.13$ at the position of galaxy 4, resulting in the $3\sigma$ upper limit of 33 $\mu$Jy at 3 GHz. Assuming radio spectral index of -0.7 (i.e., $S_\nu\propto\nu^{-0.7}$) and SFR/($M_\odot$ yr$^{-1}$)=$6.35\times10^{-22}L_\mathrm{1.4GHz}$/(W Hz$^{-1}$) \citep{2011ApJ...737...67M}, the flux upper limit translates to SFR upper limit of $\sim220$ $M_\odot$ yr$^{-1}$. Obviously, current FIR and radio data are not deep enough to constrain galaxy 4's SFR or provide informative upper limits. Future follow-up observation at submm/mm wavelengths will be critical to confirm its quiescence. If our current SFR estimate of galaxy 4 is broadly correct, this object will be particularly interesting for future studies to test galaxy quenching scenarios as discussed in the previous subsection.

\section{Summary and conclusion}\label{sec:summary}
In this work, we have analyzed the data from $\sim29$ GHz bandwidth ALMA band 3 spectral scan observation towards SMG SSA22-AzTEC26 and detected four line emitters at SNR above 5.2 in the spectral cube smoothed with 300 km s$^{-1}$ FWHM Gaussian filter, after masking the central SMG. We combine ALMA band 4/7 and multiwavelength ancillary data to derive photometric redshift and their physical properties. Our main findings are:

\par 1. Using photometric redshift, we identify two of the sources are CO$(2-1)$ at $z=1.113$ and $z=1.146$ respectively. Another is CO(3-2) at $z=2.124$ and confirmed with [CI](1-0) emission line detected in ALMA band 4 observation. The remaining source might be $z_{\mathrm{CO(2-1)}}=1.512$ or $z_{\mathrm{CO(3-2)}}=2.767$ and need to be verified by future observation.

\par
2. All three sources with CO redshift solution are massive galaxies ($M_\star>10^{10.5}M_\odot$), with SFR in the order of $\sim1-100$ $M_\odot$ yr$^{-1}$. Two of them lie within the scatter of the MS and the most massive galaxy 4 is significantly below the MS. 

\par
3. Using the CO data we estimate the molecular gas mass of $2.9-12.5\times10^{10}M_\odot$ under the assumption of normal SFG excitation and $\alpha_{\mathrm{CO}}$. We compare different gas tracers and conversion factors, finding our current choice of normal SFG excitation and $\alpha_{\mathrm{CO}}$ seems to be appropriate. The gas fraction and $t_\mathrm{dep}$ of galaxies 1 and 2 are consistent with the gas scaling relation. However, the most massive one, galaxy 4, is likely quiescent but keeps a gas ratio of $\sim27\%$, comparable with star-forming galaxies on the MS.

\par
It is difficult to explain the large gas reservoir of galaxy 4 as being leftover after quenching because of its high gas fraction and old stellar age. We speculate that the gas content is obtained in later times, via accretion or minor merger, while various quenching mechanisms might have acted to suppress star formation. Future high-resolution observations are needed to investigate stellar and gas kinematics to understand the origin of the large gas reservoir and exact quenching channels in this galaxy.

\begin{acknowledgements}
 HU acknowledges the support from JSPS KAKENHI Grant Number 20H01953. KK and RK  acknowledge the support by JSPS KAKENHI Grant Number JP17H06130. YT acknowledges the support from NAOJ ALMA Scientific Research Grant Number 2018-09B
and JSPS KAKENHI No. 17H06130, 19H01931, and 22H04939. This paper makes use of the following ALMA data: ADS/JAO.ALMA\#2019.1.01102.S and \#2021.1.01207.S. ALMA is a partnership of ESO (representing its member states), NSF (USA), and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. 
\end{acknowledgements}

\software{numpy \citep{harris2020array},
scipy \citep{2020SciPy-NMeth},
matplotlib \citep{Hunter:2007},
astropy \citep{2018AJ....156..123A}, photutils \citep{2022zndo...6825092B}, dynesty \citep{2020MNRAS.493.3132S}, SExtractor \citep{1996A&AS..117..393B}, PSFEx \citep{2011ASPC..442..435B}}

\bibliographystyle{aasjournal}


\subsubsection*{#1}}

\pagestyle{headings}
\markright{Reference sheet: \texttt{natbib}}

\usepackage{shortvrb}
\MakeShortVerb{\|}

\begin{document}
\thispagestyle{plain}

\newcommand{\textsc{Bib}\TeX}{\textsc{Bib}\TeX}
\newcommand{\texttt{#1}\def\filedate{#2}\def\fileversion{#3}}}{\texttt{#1}\def\filedate{#2}\def\fileversion{#3}}}
\begin{center}{\bfseries\Large
    Reference sheet for \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}\ usage}\\
    \large(Describing version \fileversion\ from \filedate)
\end{center}

\begin{quote}\slshape
For a more detailed description of the \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}\ package, \LaTeX\ the
source file \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}\texttt{.dtx}.
\end{quote}
\head{Overview}
The \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}\ package is a reimplementation of the \LaTeX\ |\cite| command,
to work with both author--year and numerical citations. It is compatible with
the standard bibliographic style files, such as \texttt{plain.bst}, as well as
with those for \texttt{harvard}, \texttt{apalike}, \texttt{chicago},
\texttt{astron}, \texttt{authordate}, and of course \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}.

\head{Loading}
Load with |\usepackage[|\emph{options}|]{|\texttt{#1}\def\filedate{#2}\def\fileversion{#3}}|}|. See list of
\emph{options} at the end.

\head{Replacement bibliography styles}
I provide three new \texttt{.bst} files to replace the standard \LaTeX\
numerical ones:
\begin{quote}\ttfamily
 plainnat.bst \qquad abbrvnat.bst \qquad unsrtnat.bst
\end{quote}
\head{Basic commands}
The \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}\ package has two basic citation commands, |\citet| and
|\citep| for \emph{textual} and \emph{parenthetical} citations, respectively.
There also exist the starred versions |\citet*| and |\citep*| that print
the full author list, and not just the abbreviated one.
All of these may take one or two optional arguments to add some text before
and after the citation.
\begin{quote}
\begin{tabular}{l@{\quad$\Rightarrow$\quad}l}
  |\citet{jon90}| & Jones et al. (1990)\\
  |\citet[chap.~2]{jon90}| & Jones et al. (1990, chap.~2)\\[0.5ex]
  |\citep{jon90}| & (Jones et al., 1990)\\
  |\citep[chap.~2]{jon90}| & (Jones et al., 1990, chap.~2)\\
  |\citep[see][]{jon90}| & (see Jones et al., 1990)\\
  |\citep[see][chap.~2]{jon90}| & (see Jones et al., 1990, chap.~2)\\[0.5ex]
  |\citet*{jon90}| & Jones, Baker, and Williams (1990)\\
  |\citep*{jon90}| & (Jones, Baker, and Williams, 1990)
\end{tabular}
\end{quote}
\head{Multiple citations}
Multiple citations may be made by including more than one
citation key in the |\cite| command argument.
\begin{quote}
\begin{tabular}{l@{\quad$\Rightarrow$\quad}l}
  |\citet{jon90,jam91}| & Jones et al. (1990); James et al. (1991)\\
  |\citep{jon90,jam91}| & (Jones et al., 1990; James et al. 1991)\\
  |\citep{jon90,jon91}| & (Jones et al., 1990, 1991)\\
  |\citep{jon90a,jon90b}| & (Jones et al., 1990a,b)
\end{tabular}
\end{quote}

\head{Numerical mode}
These examples are for author--year citation mode. In numerical mode, the
results are different.
\begin{quote}
\begin{tabular}{l@{\quad$\Rightarrow$\quad}l}
  |\citet{jon90}| & Jones et al. [21]\\
  |\citet[chap.~2]{jon90}| & Jones et al. [21, chap.~2]\\[0.5ex]
  |\citep{jon90}| & [21]\\
  |\citep[chap.~2]{jon90}| & [21, chap.~2]\\
  |\citep[see][]{jon90}| & [see 21]\\
  |\citep[see][chap.~2]{jon90}| & [see 21, chap.~2]\\[0.5ex]
  |\citep{jon90a,jon90b}| & [21, 32]
\end{tabular}
\end{quote}
\head{Suppressed parentheses}
As an alternative form of citation, |\citealt| is the same as |\citet| but
\emph{without parentheses}. Similarly, |\citealp| is |\citep| without
parentheses. Multiple references, notes, and the starred variants
also exist.
\begin{quote}
\begin{tabular}{l@{\quad$\Rightarrow$\quad}l}
  |\citealt{jon90}| & Jones et al.\ 1990\\
  |\citealt*{jon90}| & Jones, Baker, and Williams 1990\\
  |\citealp{jon90}| & Jones et al., 1990\\
  |\citealp*{jon90}| & Jones, Baker, and Williams, 1990\\
  |\citealp{jon90,jam91}| & Jones et al., 1990; James et al., 1991\\
  |\citealp[pg.~32]{jon90}| & Jones et al., 1990, pg.~32\\
  |\citetext{priv.\ comm.}| & (priv.\ comm.)
\end{tabular}
\end{quote}
The |\citetext| command
allows arbitrary text to be placed in the current citation parentheses.
This may be used in combination with |\citealp|.
\head{Partial citations}
In author--year schemes, it is sometimes desirable to be able to refer to
the authors without the year, or vice versa. This is provided with the
extra commands
\begin{quote}
\begin{tabular}{l@{\quad$\Rightarrow$\quad}l}
  |\citeauthor{jon90}| & Jones et al.\\
  |\citeauthor*{jon90}| & Jones, Baker, and Williams\\
  |\citeyear{jon90}|   & 1990\\
  |\citeyearpar{jon90}| & (1990)
\end{tabular}
\end{quote}
\head{Forcing upper cased names}
If the first author's name contains a \textsl{von} part, such as ``della
Robbia'', then |\citet{dRob98}| produces ``della Robbia (1998)'', even at the
beginning of a sentence. One can force the first letter to be in upper case
with the command |\Citet| instead. Other upper case commands also exist.
\begin{quote}
\begin{tabular}{rl@{\quad$\Rightarrow$\quad}l}
  when & |\citet{dRob98}| & della Robbia (1998) \\
  then & |\Citet{dRob98}| & Della Robbia (1998) \\
     &   |\Citep{dRob98}| & (Della Robbia, 1998) \\
     &   |\Citealt{dRob98}| & Della Robbia 1998 \\
     &   |\Citealp{dRob98}| & Della Robbia, 1998 \\
     &   |\Citeauthor{dRob98}| & Della Robbia
\end{tabular}
\end{quote}
These commands also exist in starred versions for full author names.

\head{Citation aliasing}
Sometimes one wants to refer to a reference with a special designation,
rather than by the authors, i.e. as Paper~I, Paper~II. Such aliases can be
defined and used, textual and/or parenthetical with:
\begin{quote}
\begin{tabular}{lcl}
  |\defcitealias{jon90}{Paper~I}|\\
  |\citetalias{jon90}| & $\Rightarrow$ & Paper~I\\
  |\citepalias{jon90}| & $\Rightarrow$ & (Paper~I)
\end{tabular}
\end{quote}
These citation commands function much like |\citet| and |\citep|: they may
take multiple keys in the argument, may contain notes, and are marked as
hyperlinks.
\head{Selecting citation style and punctuation}
Use the command |\bibpunct| with one optional and 6 mandatory arguments:
\begin{enumerate}
\item the opening bracket symbol, default = (
\item the closing bracket symbol, default = )
\item the punctuation between multiple citations, default = ;
\item the letter `n' for numerical style, or `s' for numerical superscript
      style, any other letter for
      author--year, default = author--year;
\item the punctuation that comes between the author names and the year
\item the punctuation that comes between years or numbers when common author
      lists are suppressed (default = ,);
\end{enumerate}

The optional argument is the character preceding a post-note, default is a
comma plus space. In redefining this character, one must include a space if
one is wanted.

Example~1, |\bibpunct{[}{]}{,}{a}{}{;}| changes the output of
\begin{quote}
  |\citep{jon90,jon91,jam92}|
\end{quote}
into [Jones et al. 1990; 1991, James et al. 1992].

Example~2, |\bibpunct[; ]{(}{)}{,}{a}{}{;}| changes the output of
\begin{quote}
  |\citep[and references therein]{jon90}|
\end{quote}
into (Jones et al. 1990; and references therein).

\head{Other formatting options}
Redefine |\bibsection| to the desired sectioning command for introducing
the list of references. This is normally |\section*| or |\chapter*|.

Define |\bibpreamble| to be any text that is to be printed after the heading but
before the actual list of references.

Define |\bibfont| to be a font declaration, e.g.\ |\small| to apply to
the list of references.

Define |\citenumfont| to be a font declaration or command like |\itshape|
or |\textit|.

Redefine |\bibnumfmt| as a command with an argument to format the numbers in
the list of references. The default definition is |[#1]|.

The indentation after the first line of each reference is given by
|\bibhang|; change this with the |\setlength| command.

The vertical spacing between references is set by |\bibsep|; change this with
the |\setlength| command.

\head{Automatic indexing of citations}
If one wishes to have the citations entered in the \texttt{.idx} indexing
file, it is only necessary to issue |\citeindextrue| at any point in the
document. All following |\cite| commands, of all variations, then insert
the corresponding entry to that file. With |\citeindexfalse|, these
entries will no longer be made.

\head{Use with \texttt{chapterbib} package}

The \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}\ package is compatible with the \texttt{chapterbib} package
which makes it possible to have several bibliographies in one document.

The package makes use of the |\include| command, and each |\include|d file
has its own bibliography.

The order in which the \texttt{chapterbib} and \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}\ packages are loaded
is unimportant.

The \texttt{chapterbib} package provides an option \texttt{sectionbib}
that puts the bibliography in a |\section*| instead of |\chapter*|,
something that makes sense if there is a bibliography in each chapter.
This option will not work when \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}\ is also loaded; instead, add
the option to \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}.

Every |\include|d file must contain its own
|\bibliography| command where the bibliography is to appear. The database
files listed as arguments to this command can be different in each file,
of course. However, what is not so obvious, is that each file must also
contain a |\bibliographystyle| command, \emph{preferably with the same
style  argument}.
\head{Sorting and compressing citations}
Do not use the \texttt{cite} package with \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}; rather use one of the
options \texttt{sort} or \texttt{sort\&compress}.

These also work with author--year citations, making multiple citations appear
in their order in the reference list.

\head{Long author list on first citation}
Use option \texttt{longnamesfirst} to have first citation automatically give
the full list of authors.

Suppress this for certain citations with |\shortcites{|\emph{key-list}|}|,
given before the first citation.

\head{Local configuration}
Any local recoding or definitions can be put in \texttt{#1}\def\filedate{#2}\def\fileversion{#3}}\texttt{.cfg} which
is read in after the main package file.

\head{Options that can be added to \texttt{\char`\\ usepackage}}
\begin{description}
\item[\ttfamily round] (default) for round parentheses;
\item[\ttfamily square] for square brackets;
\item[\ttfamily curly] for curly braces;
\item[\ttfamily angle] for angle brackets;
\item[\ttfamily colon] (default) to separate multiple citations with
     colons;
\item[\ttfamily comma] to use commas as separaters;
\item[\ttfamily authoryear] (default) for author--year citations;
\item[\ttfamily numbers] for numerical citations;
\item[\ttfamily super] for superscripted numerical citations, as in
     \textsl{Nature};
\item[\ttfamily sort] orders multiple citations into the sequence in
     which they appear in the list of references;
\item[\ttfamily sort\&compress] as \texttt{sort} but in addition multiple
     numerical citations are compressed if possible (as 3--6, 15);
\item[\ttfamily longnamesfirst] makes the first citation of any reference
     the equivalent of the starred variant (full author list) and subsequent
     citations normal (abbreviated list);
\item[\ttfamily sectionbib] redefines |\thebibliography| to issue
     |\section*| instead of |\chapter*|; valid only for classes with a
     |\chapter| command; to be used with the \texttt{chapterbib} package;
\item[\ttfamily nonamebreak] keeps all the authors' names in a citation on
     one line; causes overfull hboxes but helps with some \texttt{hyperref}
     problems.
\end{description}
\end{document}
