# Deep Observations of CO and Free-Free Emission in Ultraluminous Infrared QSO
IRAS F07599+6508

Qing-Hua Tan Purple Mountain Observatory & Key Laboratory for Radio Astronomy,
Chinese Academy of Sciences, 10 Yuanhua Road, Nanjing 210023, People’s
Republic of China<EMAIL_ADDRESS>Yu Gao Department of Astronomy, Xiamen
University, 422 Siming South Road, Xiamen 361005, People’s Republic of China;
<EMAIL_ADDRESS>Purple Mountain Observatory & Key Laboratory for Radio
Astronomy, Chinese Academy of Sciences, 10 Yuanhua Road, Nanjing 210023,
People’s Republic of China<EMAIL_ADDRESS>Emanuele Daddi CEA, Irfu, DAp,
AIM, Universitè Paris-Saclay, Universitè de Paris, CNRS, F-91191 Gif-sur-
Yvette, France  Xiao-Yang Xia  Tianjin Astrophysics Center, Tianjin Normal
University, Tianjin 300387, People’s Republic of China  Cai-Na Hao Tianjin
Astrophysics Center, Tianjin Normal University, Tianjin 300387, People’s
Republic of China  Alain Omont Sorbonne Universit$\acute{e}$, UPMC
Universit$\acute{e}$ Paris 6 and CNRS, UMR 7095, Institut d’Astrophysique de
Paris, 98bis boulevard Arago, 75014 Paris, France  Kotaro Kohno Institute of
Astronomy, Graduate School of Science, The University of Tokyo, Osawa, Mitaka,
Tokyo 181-0015, Japan  Research Center for the Early Universe, Graduate School
of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033,
Japan

###### Abstract

Infrared quasi-stellar objects (IR QSOs) are a rare subpopulation selected
from ultraluminous infrared galaxies (ULIRGs) and have been regarded as
promising candidates of ULIRG-to-optical QSO transition objects. Here we
present NOEMA observations of the CO(1$-$0) line and 3 mm continuum emission
in an IR QSO IRAS F07599+6508 at $z=0.1486$, which has many properties in
common with Mrk 231. The CO emission is found to be resolved with a major axis
of $\sim$6.1 kpc that is larger than the size of $\sim$4.0 kpc derived for 3
mm continuum. We identify two faint CO features located at a projected
distance of $\sim$11.4 and 19.1 kpc from the galaxy nucleus, respectively,
both of which are found to have counterparts in the optical and radio bands
and may have a merger origin. A systematic velocity gradient is found in the
CO main component, suggesting that the bulk of molecular gas is likely
rotationally supported. Based on the radio-to-millimeter spectral energy
distribution and IR data, we estimate that about 30$\%$ of the flux at 3 mm
arises from free-free emission and infer a free-free-derived star formation
rate of 77 $M_{\odot}\ {\rm yr^{-1}}$, close to the IR estimate corrected for
the AGN contribution. We find a high-velocity CO emission feature at the
velocity range of about $-$1300 to $-$2000 km s-1. Additional deep CO
observations are needed to confirm the presence of a possible very high-
velocity CO extension of the OH outflow in this IR QSO.

galaxies: active - galaxies: evolution - galaxies: ISM - galaxies: starburst -
ISM: molecules

††facilities: IRAM:NOEMA, VLA, Sloan, HST, PS1††software: GILDAS (Guilloteau &
Lucas, 2000), CASA (McMullin et al., 2007), AIPS (van Moorsel et al., 1996)

##  1 Introduction

Mergers of gas-rich galaxies may trigger both intense star formation and
nuclear activity, leading to a violent transition from an obscured accretion
and starburst stage to an optically bright quasi-stellar object (QSO) and
eventually an elliptical galaxy (e.g., Sanders et al., 1988a, b; Barnes &
Hernquist, 1992; Di Matteo et al., 2005; Hopkins et al., 2006). In this
merger-driven evolutionary scenario, the active galactic nucleus (AGN) was
thought to be energetically dominant but obscured in the ultraluminous
infrared galaxy (ULIRG; $L_{\rm IR}\geqslant 10^{12}\ L_{\odot}$) phase, where
the most luminous sources that reach a QSO-like luminosity represent the rare
ULIRG-to-QSO transition objects. These systems are very likely experiencing
the so-called feedback or blowout phase that expels the gas of the host galaxy
in the form of outflowing winds, which have been invoked as one of the main
drivers for regulating the growth of the stellar mass and black hole (BH)
accretion and producing the tight BH-bulge mass relation (see Silk & Rees,
1998; Hopkins et al., 2008; Fabian, 2012; Kormendy & Ho, 2013; King & Pounds,
2015, and references therein), and therefore represent a key phase in the
lifetime of a quasar.

The large gas masses discovered in local ULIRGs via CO observations have
established the presence of abundant fuel for the intense star formation
activity (e.g., Solomon et al., 1997). Along with the clear evidence that many
(U)LIRGs are strong dynamical interactions and mergers (Sanders & Mirabel,
1996), all the properties make these systems unique local examples of dust-
enshrouded galaxies at high-$z$. Xia et al. (2012) have investigated for the
first time the CO survey in 19 ultraluminous infrared type 1 QSOs, selected
from the local ULIRGs and referred to as IR QSOs (see Zheng et al., 2002).
These IR-selected type 1 QSOs (IR QSOs) were originally compiled from the
ULIRGs in the QDOT redshift survey (Lawrence et al., 1999), the 1 Jy ULIRG
survey (Kim & Sanders, 1998), and the cross-correlation of the IRAS Point
Source Catalog with the ROSAT All-sky Survey Catalog. Studies of luminosity
function for large samples of IRAS galaxies reveal that the space densities of
ULIRGs and QSOs are comparable in the local universe (Soifer et al., 1987),
while the fraction of IR QSOs is less than 10% of ULIRGs (Lawrence et al.,
1999), indicating that the IR QSO phase may last only a few times 107 yr if
the space density of objects is simply related to the timescale of different
phases (Hao et al., 2005).

The CO observations of IR QSOs show that the molecular gas masses are a few
times $10^{9}-10^{10}\ M_{\odot}$ (Xia et al., 2012), similar to that found in
local ULIRGs but significantly higher than that of optically selected PG QSOs
(e.g., Scoville et al., 2003; Evans et al., 2001, 2006; Shangguan et al.,
2020). This supports the scenario that IR QSOs are rich of fuel to feed both
the star formation and the accretion of AGN. Similar conclusion has also been
reached recently for type 2 QSOs, of which a sample of QSO hosts are found to
be gas-rich based on the CO observations (Jarvis et al., 2020). Together with
the far-IR excess observed in IR QSOs (Hao et al., 2005), and the multi-
wavelengths properties (i.e., large blueshifts in optical emission lines, warm
IR color, intermediate slope of mid-IR spectral energy distribution (SED)
between ULIRGs and PG QSOs, and high accretion rates; Zheng et al., 2002; Cao
et al., 2008), all of these properties imply that IR QSOs are likely the
objects caught in the short-lived transition phase between ULIRG and QSO
stages (Canalizo & Stockton, 2001).

Compared to the local (U)LIRGs that a large number have been spatially
resolved in molecular gas emission (e.g., Downes & Solomon, 1998; Wilson et
al., 2008; Ueda et al., 2014; Xu et al., 2014, 2015; Imanishi et al., 2019),
except Mrk 231, the nearest IR QSO at $z=0.042$ that has been well-studied in
both multi-wavelength and multi-scale (e.g., Bryant & Scoville, 1996; Taylor
et al., 1999; Fischer et al., 2010; Rupke & Veilleux, 2011; Feruglio et al.,
2015; Cicone et al., 2020), until recently the CO emission have been spatially
resolved on sub-kiloparsec (sub-kpc) scale in a handful of IR QSOs with the
ALMA (Tan et al., 2019). It was found that the host galaxies of IR QSOs
resolved with ALMA have a diverse variety of CO morphology and kinematics
(i.e., rotating disk, disturbed structure, and multiple interacting sources),
indicative of more complicated evolutionary stages from merging (U)LIRGs to
QSOs. This is somewhat unexpected since the objects in IR QSO phase are
typically characterized by the final coalescence of the galaxies, according to
the models of merger-driven evolutionary scenario (Barnes & Hernquist, 1992;
Hopkins et al., 2008). However, such a classification suffers from small
number statistics, spatially resolved CO observations in more IR-luminous QSOs
are required for a meaningful statistical study.

In this paper, we present CO(1$-$0) imaging of IRAS F07599+6508 obtained with
the IRAM NOrthern Extended Millimeter Array (NOEMA) interferometer. IRAS
F07599+6508 is an IR QSO at $z=0.149$ with a total IR luminosity of $3.5\times
10^{12}\ L_{\odot}$ and AGN bolometric luminosity of $L_{\rm
AGN,bol}=3.9\times 10^{12}\ L_{\odot}$ (Hao et al., 2005). It is an active
galaxy of particular interest because of its highly similarity to Mrk 231, the
closest and most prominent template of local IR QSOs, but at a higher redshift
and lack of detailed studies on the cold molecular gas of the host galaxy,
especially the spatial distribution and kinematics. Similar to Mrk 231, IRAS
F07599+6508 is a low-ionization broad absorption-line (BAL) QSO with strong
Feii emission (e.g., Kwan et al., 1995; Véron-Cetty et al., 2006), which is a
very rare subclass of BAL QSOs that may represent young, heavily enshrouded
AGNs, where the cocoon of gas and dust has a large covering factor (Boroson &
Meyers, 1992; Zheng et al., 2002).

The IRAM 30m CO observations show massive amounts of molecular gas residing in
IRAS F07599+6508 and a broad FWHM line width of $380\pm 14$ km s-1 and $490\pm
40$ km s-1 for the CO(1$-$0) and CO(2$-$1) line, respectively (Xia et al.,
2012). Based on six independent methods, Veilleux et al. (2009) estimate an
average fraction of AGN contribution to the bolometric luminosity of
$\alpha_{\rm AGN}=0.75$, consistent with other IR QSOs (e.g., $\sim 70\%$ for
Mrk 231) where the $\alpha_{\rm AGN}$ is found to be typically higher than
that ($\sim 35\%-40\%$) of local cooler ULIRGs. The weakness of both soft
($L_{\rm 0.1-2.4\ keV}=6.9\times 10^{42}\ {\rm erg\ s^{-1}}$) and hard X-ray
emission ($L_{\rm 2-10\ keV}=1.3\times 10^{42}\ {\rm erg\ s^{-1}}$) seen in
IRAS F07599+6508 (e.g., Xia et al., 2001; Brightman & Nandra, 2011; Luo et
al., 2014; La Caria et al., 2019), i.e., the X-ray emission is not at the
level expected from its optical/UV emission compared to typical QSOs, imply
that its central AGN is very likely highly obscured.

The existence of multi-phase outflow winds in IRAS F07599+6508 has been
demonstrated by the detection of large blueshifts ($\sim$2000 km s-1) of
optical emission lines, Na i D interstellar absorption line, and OH 119$\mu$m
P-Cygni/blue-shifted absorption profile, closely resemble to Mrk 231 (Boroson
& Meyers, 1992; Zheng et al., 2002; Spoon et al., 2013; Veilleux et al., 2013;
Rupke et al., 2017). The sensitive CO(1$-$0) mapping observations will help
constrain the distribution and kinematics of the molecular gas reservoir in
IRAS F07599+6508 and offer the possibility of detecting the molecular outflow
based on the presence of broad wings in the CO line emission. Tremendous
amount of observing time were devoted to Mrk 231 on NOEMA revealing the
outstanding molecular outflows including in dense gas tracers (Aalto et al.,
2012), and likely the tentative detection of large-scale molecular oxygen
emission (Wang et al., 2020).

Alike to Mrk 231, IRAS F07599+6508 has typical features expected for a system
transiting from the heavily obscured starburst phase to the unobscured QSO
phase, and could also be regarded as a representative of IR QSOs which have
low number density in the local universe but are expected to be common at high
redshift (Veilleux et al., 1999; Le Floc’h et al., 2005), and thus provides
nearby laboratory to probe the ULIRG-to-optical QSO transition phase that may
hold the key to understanding the galaxy evolution in a merger-driven
scenario.

The paper is organized as follows. In Section 2, we describe the NOEMA
observations and data reduction process. Section 3 presents the results from
the analysis of 3 mm continuum and CO(1$-$0) line data, including morphology,
spatially decomposed CO components, kinematics, and an analysis of archival
radio continuum data from VLA. In section 4, we discuss the nature of 3 mm
continuum emission, the origin of the off-center CO components, and make a
comparison of the molecular gas properties between IRAS F07599+6508 and Mrk
231. We summarize our main findings in Section 5. Throughout this work, we
assume a standard $\Lambda$CDM cosmology with $H_{\rm 0}$ = 70 km s-1 Mpc-1,
$\Omega_{\rm m}$ = 0.3, and $\Omega_{\Lambda}$ = 0.7.

##  2 Observations and data reduction

###  2.1 NOEMA CO(1-0) Observations

We have observed the CO(1$-$0) line and the continuum emission in IRAS
F07599+6508 with the IRAM NOEMA interferometer between May and July 2018. The
observations were obtained with a total on-source time of 12.7 h in the
compact array D-configuration (with eight antennas covering baselines of 15.2
- 175.9 m). We used the PolyFix correlator in 3 mm Band 1, which provides
2$\times$8 GHz of instantaneous dual-polarisation bandwidth with a spectral
resolution of 2 MHz (corresponding to 6.0 km s-1at the redshifted CO(1$-$0)
frequency of 100.353 GHz). The tuning frequency was set to 98.2 GHz with
frequency ranged from 94.2 to 102.2 GHz in the Upper SideBand (USB) and 78.8
to 86.8 GHz in the Lower SideBand (LSB). In this paper, we analyze the CO(1-0)
line and the continuum data in the USB; the rest of detected lines (i.e.,
CN(1-0) lines at 113.191 and 113.491 GHz, CO(1$-$0) isotopologue lines, and
the 109.174 GHz HC3N (12$-$11) line) covered in the USB and the LSB line data
will be presented in a future paper.

The data calibration was performed using the
GILDAS111http://iram.fr/IRAMFR/GILDAS; for more information, see Guilloteau &
Lucas (2000) package CLIC with help from the staff in Grenoble. The typical
uncertainty of flux calibration in the 3 mm band is about 10%. We binned the
$uv$ tables by 20 km s-1and redefined the rest frequency to the sky frequency
of CO(1$-$0) transition. We then used the software GILDAS/MAPPING for image
cleaning and analysis in the $uv$-plane. We also exported the data to a
CASA222http://casa.nrao.edu; see McMullin et al. (2007)-compatible format for
a comparison analysis. The continuum $uv$ table was created from the
calibrated visibilities by averaging the line-free channels using the MAPPING
task uv_continuum with line emission channels filtered (uv_filter task). The
continuum emission was subtracted from the line $uv$-data by fitting a first-
order polynomial baseline for each visibility (uv_baseline task). The USB and
LSB continuum $uv$-data were combined using the uv_merge task to form a single
output table centered at 90.5 GHz. The synthesized beam size using natural
weighting is 5.″67$\times$3.″81 and the NOEMA primary beam at the sky
frequency of CO(1$-$0) is 50.″2$\times$50.″2. We obtain a sensitivity of
$\sim$0.44 mJy beam-1 in the spectral velocity resolution of 20 km s-1 and
$\sim$10 $\mu$Jy beam-1 for the continuum image.

###  2.2 Archival VLA Radio Data

In this work we also make use of ancillary radio continuum data observed with
the Very Large Array (VLA) for a comprehensive analysis. IRAS F07599+6508 has
been observed at 1.49, 4.86, 8.44, 14.94, and 22.46 GHz in CD configuration
(project ID: AB0783). Higher resolution observations in AB configuration
(project ID: AH0406) were available at 1.49, 4.86, and 8.44 GH, as well as in
A configuration at 14.94 GHz (project ID: AN0104). We obtained the raw data
sets from the VLA
archive333https://science.nrao.edu/facilities/vla/archive/index. The data
flagging and reduction were performed by using the AIPS package (van Moorsel
et al., 1996). The full image at each frequency and array was processed
individually with the AIPS CLEAN task imagr and then transferred to FITS
format for further analysis with the CASA package. The synthesized beam using
Briggs weighting with a robustness parameter of 0.5 ranges from
65.″2$\times$28.″3 at 1.49 GHz to 1.″67$\times$1.″30 at 43.34 GHz for
observations in CD configuration, and from 5.″33$\times$2.″08 at 1.49 GHz to
0.″18$\times$0.″13 at 14.94 GHz for higher resolution observations in AB and A
configurations.

##  3 Results and analysis

###  3.1 3 mm Continuum Emission

Figure 1: NOEMA cleaned map of 3 mm continuum emission (left), the residual
map of 3 mm emission after subtracting a circular Gaussian model (middle), and
the corresponding $uv$-plot showing the visibility amplitude as a function of
the $uv$-distance (right) for IRAS F07599+6508. For 3 mm continuum map,
contours start at $\pm 2\sigma$ (1$\sigma$ rms noise is 10 $\mu$Jy beam-1) and
increase by a factor of 1.5, with positive (negative) contours shown as solid
(dotted) lines, while for the residual map, contours start at $\pm 2\sigma$
and increase in steps of 1$\sigma$. The synthesized beam of 6.″25$\times$4.″22
with position angle 106.2∘ is shown in the bottom left corner. The plus sign
indicates the centroid of the 3 mm continuum emission derived from a circular
Gaussian fit to the $uv$-data. In the right panel, the black circles with
error bars shown in the $uv$-plot are the data binned in $uv$-radius steps of
8 k$\lambda$. Error bars show the statistical photon noise on the average
amplitudein each bin. The red curve indicates the best-fitting model with a
circular Gaussian function to the $uv$-data.

The continuum emission at 90.5 GHz (rest-frame 103.9 GHz; $\sim$3 mm) is
clearly detected towards the nucleus of IRAS F07599+6508. Figure 1 shows the 3
mm continuum map and the corresponding $uv$-plot that shows the amplitude of
the visibilities as a function of the $uv$-distance. It can be seen from the
$uv$-data that the 3 mm continuum source is resolved in our observations, as
an unresolved source would give a flat distribution of visibilities.

Assuming that the surface brightness distribution of the continuum follow a
symmetric Gaussian profile, we performed a fit to the visibilities in GILDAS
using task uv_fit and found that the continuum source is best fitted by a two-
dimensional (2D) circular Gaussian model with a deconvolved full width at half
maximum (FWHM) size of 1.″53$\pm$0.″19 ($\sim 4.0\pm 0.5$ kpc) and a flux
density of 560$\pm$20 $\mu$Jy (see Table 1). A 2D elliptical Gaussian model
fit to the continuum visibilities returns a deconvolved FWHM size of
(1.″54$\pm$0.″22)$\times$(1.″50$\pm$0.″39) and a flux density of 560$\pm$20
$\mu$Jy, which are in good agreement with those derived from fit with a
circular Gaussian model, but with slightly larger uncertainties on the size
measurement. In addition, the position angle ($0^{\circ}\pm 293^{\circ}$)
derived from elliptical Gaussian fit is rather uncertain due to the large
uncertainty.

In the middle panel of Figure 1, we show the residual after removing the fit
with a circular Gaussian model. No significant peaks or negative regions with
signal-to-noise ratio (SNR) $\gtrsim 3$ are found in the residual map. The red
curve in $uv$-plot (see the right panel of Figure 1) shows the result of
circular Gaussian model fit. Moreover, it can be seen from the right panel of
Figure 1, most data points are consistent with the model fit (see the red
curve) with a circular Gaussian function within 1$\sigma$. All of these
suggest that a single circular Gaussian model may provide a good description
to the 3 mm continuum emission based on the current data.

As an independent check on the results from our analysis of $uv$-data, we also
measured the size and flux of the continuum emission by fitting the data in
the image plane using the imfit task in CASA. A 2D Gaussian fit to the 3 mm
continuum map of Figure 1 finds an deconvolved FWHM size of
(1.″84$\pm$0.″45)$\times$(1.″40$\pm$0.″49) and an integrated flux of
550$\pm$20 $\mu$Jy, consistent with the estimates from the uv_fit within
1$\sigma$ uncertainties.

Figure 2: NOEMA continuum-subtracted spectrum (left) and cleaned map (right)
of CO(1$-$0) emission of IRAS F07599+6508. The CO(1$-$0) spectrum is extracted
from the central beam and the spectrum is binned to a channel width of
$\sim$20 km s-1. The red dotted line represents a Gaussian fit to the
spectrum. The zero velocity corresponds to the redshift of $z_{\rm CO}=0.1486$
derived from a Gaussian fit to the spectrum. The map of CO(1$-$0) emission is
integrated in the [$-$630, +450] km s-1 velocity range. Contours start at $\pm
2\sigma$ (1$\sigma$ rms noise is 0.092 Jy beam-1 km s-1) and increase by a
factor of 1.5, with positive (negative) contours shown as solid (dotted)
lines. The synthesized beam (5.″61$\times$3.″77) is shown in the bottom left
corner of the map. The plus sign indicates the peak of 3 mm continuum
emission.

###  3.2 CO(1$-$0) Line Emission

####  3.2.1 CO morphology and spatial decomposition

Figure 2 shows the continuum-subtracted CO(1$-$0) spectrum and the integrated
intensity map of IRAS F07599+6508. The CO spectrum is extracted by fitting a
point-source model to the $uv$-table using the uv_fit task. By fitting a
Gaussian profile to the spectrum, we measure a redshift of $z_{\rm
CO}=0.1486$, being used to set the systematic velocity. The velocity width
(FWHM) derived from the Gaussian fit to the spectrum is 430$\pm$10 km s-1. We
note that the CO(1$-$0) line profile is asymmetrical and a single Gaussian
function does not provide a good fit to the line peak and wings. A double-peak
profile is detected at the line center albeit not very prominent, which is
often taken to indicate orbital motion in a merger or a rotating disk (e.g.,
Downes & Solomon, 1998; Greve et al., 2005; Daddi et al., 2010). In addition,
the line profile also shows a slightly weaker peak at velocity $\sim 100-150$
km s-1, indicating that the motion of molecular gas in the host galaxy is
complex.

The irregular morphology of the CO emission with extended structure is shown
by the map in the right panel of Figure 2, obtained by collapsing the channels
in the [$-$630, +450] km s-1 velocity range where the flux is positive in each
channel when the CO spectrum is rebinned to a channel width of $\sim 60$ km
s-1. The CO emission appears to be resolved with an extension to the northeast
and to the south of the nucleus, respectively. To characterize the
distribution of CO emission, we fitted an elliptical Gaussian model to the
$uv$-data with channels averaged in the [$-$630, +450] km s-1 velocity range
using the task uv_fit. The best fit gives a source size FWHM of $(2.\arcsec
33\pm 0.\arcsec 11)\times(1.\arcsec 91\pm 0.\arcsec 07)$. The peak of CO
emission is found to be co-spatial with the 3 mm continuum emission.

In the residual map after subtraction of an elliptical Gaussian component, we
detect two off-center CO subcomponents at above 5$\sigma$, which we denote as
northeast (NE) and south (S), and we refer to the main component fitted with
an elliptical Gaussian in the galaxy center as F07599C (see Figure 3). To take
into account the off-center components, we remake the fitting in the
$uv$-plane using combinations of functions by allowing all parameters to be
free in the modeling, and the details of the model fitting are described in
Appendix A. The best fit shows that the CO emission can be spatially
decomposed into three components, an elliptical Gaussian source (F07599C) and
two point sources (F07599S and F07599NE).

For the main component F07599C, the model fit returns a source size FWHM of
$(2.\arcsec 35\pm 0.\arcsec 11)\times(1.\arcsec 66\pm 0.\arcsec 09)$,
corresponding to 6.1 by 4.3 kpc, which is about two times of the median source
size ($\sim$3.2 kpc) measured for the local IR QSOs with well-resolved CO
morphology (Tan et al., 2019). The three components identified in the CO
emission and their parameters are listed in Table 1.

Figure 3: Residual maps of CO(1$-$0) emission after subtracting a single
elliptical Gaussian model (left) and a combination of an ellipical Gaussian
model and two point sources models (right). The black plus sign marks the peak
of 3 mm continuum emission, while the two red crosses indicate the peak
positions of the two point-like CO feature. Contours start at $\pm 2\sigma$
and increase in steps of 1$\sigma$ for both maps. The two point-like
subcomponents are labeled as F07599S and F07599NE, respectively.  Table 1:
Measured properties of CO(1$-$0) line and 3 mm continuum emission of IRAS
F07599+6508

Emission | Component | R.A. | Decl. | Size FWHM | FWHMCO(1-0) | $v_{\rm offset}$ | Integrated Flux
---|---|---|---|---|---|---|---
|  | (J2000) | (J2000) | (arcsec) | (km s-1) | (km s-1) |
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8)
CO(1$-$0) | F07599C | 08:04:30.468$\pm$0.004 | 64:59:52.72$\pm$0.02 | (2.35$\pm$0.11)$\times$(1.66$\pm$0.09) | 430$\pm$10 | 0 | 11.24$\pm$0.09 Jy km s-1
| F07599S | 08:04:30.41$\pm$0.03 | 64:59:48.36$\pm$0.12 |  $\cdots$ | 100$\pm$10 (narrow) | 100$\pm$10 | 1.12$\pm$0.07 Jy km s-1
|  |  |  |  | 500$\pm$80 (broad) |  |
| F07599NE | 08:04:31.56$\pm$0.07 | 64:59:55.2$\pm$0.3 |  $\cdots$ | 390$\pm$70 | 40$\pm$30 | 0.45$\pm$0.07 Jy km s-1
3 mm | total | 08:04:30.432$\pm$0.010 | 64:59:52.71$\pm$0.04 | 1.53$\pm$0.19 |  |  | 560$\pm$20 $\mu$Jy

Note. — Column (1): CO(1$-$0) line and 3 mm continuum emission. Column (2):
components of CO(1$-$0) emission from source spatial decomposition in the
$uv$-plane. Column (3) and (4): source position derived from model fit to the
$uv$-data. Column (5): FWHM source size derived from model fitting to the
visibilities using the GILDAS task uv_fit. Dots in the column indicate an
unresolved component. The physical scale of 1″corresponds to $\sim$2.6 kpc.
Column (6): FWHM CO line width derived from a Gaussian profile fit to the
spectrum (except for F07599S, where the CO spectrum is fitted simultaneously
with two Gaussians consist of one narrow and one broad components centered at
the same velocity). Column(7): CO(1$-$0) line centroid velocity offset from
the systematic velocity. The centroid velocity is derived from a Gaussian fit
to the spectrum. Column (8): CO(1$-$0) velocity-integrated intensity and 3 mm
continuum emission flux. The CO integrated fluxes are measured from the
spectra shown in Figure 4 by integrating the intensity over the line-emitting
region in each channel. The continuum flux is measured from a circular
Gaussian fit to the $uv$-data. Errors in fluxes quoted in this table are
statistical and do not include the absolute flux calibration uncertainty.

Figure 4: Top: CO(1$-$0) spectra of the three components extracted from GILDAS
task uv_fit using an elliptical Gaussian and two point-source functions for
F07599C (orange), F07599S (purple) and F07599NE (green), respectively. The CO
spectrum of F07599C is divided by 5 for comparison purpose. All spectra are
smoothed to a velocity resolution of $\sim$20 km s-1. Bottom-left: zoomed-in
views of CO(1$-$0) spectra for F07599S and F07599NE and the spectra are re-
binned in channels of $\sim$60 km s-1. The black lines show the Gaussian fits
to the spectra. For F07599S, we modeled the line profile with two Gaussians
(one narrow and one broad), while for F07599NE the spectrum is fitted by using
a single Gaussian. Bottom-right: the corresponding continuum-subtracted
CO(1$-$0) maps (colored contours) of the three components overlaid on the
CO(1$-$0) map of the whole IRAS F07599+6508 (greyscale; same as the map in the
right panel of Figure 2). The colored contours show the CO(1$-$0) emission
integrated in the [-330, +330] km s-1 (orange), [-270, +330] km s-1 (purple),
and [-330, +270] km s-1 (green) velocity range, respectively, while the black
contour indicates the global CO(1$-$0) emission detected at a 3$\sigma$
significance in the intensity map integrated in the [$-$630, +450] km s-1
velocity range. The colored contours start at $2\sigma$ and increase by a
factor of 1.5. The 1$\sigma$ rms noise is 0.075, 0.095, and 0.082 Jy beam-1 km
s-1 for F07599C, F07599S, and F07599NE, respectively. The CO map of F07599C
was derived from a cleaned data cube, while the maps of F07599S and F07599NE
were not cleaned given the relatively low SNR.

####  3.2.2 Properties of the decomposed CO components

Figure 4 shows the CO(1$-$0) spectra of the three components and the
corresponding velocity-integrated CO maps. The CO spectra of the three
components are extracted simultaneously from the line cube in the $uv$-plane
with fitting models identical to the best-fit models summarized in Table 1
(see Section 3.2.1). The spectrum of F07599C is obtained with $uv$-modeling
fit using an elliptical Gaussian function by fixing the parameters of centroid
and major/minor axes, while for the off-center components, F07599S and
F07599NE, which are offset from the nucleus of main component with a projected
distance of 4.″4 ($\sim$ 11.4 kpc) and 7.″4 ($\sim$ 19.1 kpc), respectively,
the spectra are extracted by a point source model with centroids fixed.

A single Gaussian fit to the line profile of F07599C and F07599NE shown in
Figure 4 yields a peak flux density of 26.0$\pm$0.3 mJy and 1.0$\pm$0.2 mJy,
and an FWHM of 430$\pm$10 km s-1 and 390$\pm$70 km s-1, respectively. F07599NE
has a velocity offset of 40$\pm$30 km s-1, consistent within the uncertainties
with the systematic velocity. For F07599S, the spectrum is best fitted by two
nested Gaussians (with a shared centre), e.g., a narrow Gaussian to fit the
central core of the line, and a broad Gaussian to fit the broad emission. The
best-fit gives a peak flux density of 3.9$\pm$0.4 mJy and 1.3$\pm$0.3 mJy, and
an FWHM of 100$\pm$10 km s-1 and 500$\pm$80 km s-1 with velocity centroid
redshifted at a velocity of about 100$\pm$10 km s-1 for the narrow and broad
emission, respectively. The results of Gaussian fits derived for the three
components are reported in Table 1.

Similar to the CO spectrum extracted from the central beam, the CO line
profile of F07599C shows characteristic double peak at the line center and a
weaker peak at the velocity $\sim 100-150$ km s-1, which is found to be
coincident with the CO line centroid velocity of F07599S (see Figure 4, top
panel). The velocity coincidence is likely to suggest that a significant
portion of molecular gas are detected in the overlap region between F07599C
and F07599S, where the perturbation of molecular gas make an imprint on the CO
line profile. Higher spatial resolution observations are needed to resolve the
emission and clarify the nature of the source.

The CO morphology of each component is shown by the map in the bottom-right
panel of Figure 4, obtained by averaging the velocity range of [$-$330, +330]
km s-1, [$-$270, +330] km s-1, and [$-$330, +270] km s-1, where the emission
is detected in consecutive channels at above 1$\sigma$, for F07599C, F07599S,
and F07599NE, respectively. The velocity-integrated CO emission of each
component was imaged separately by collapsing the velocity channels over the
line-emitting region. It can be seen from Figure 4 that the superimposed image
of the three components agrees well with the global CO map that is integrated
in the [$-$630, +450] km s-1 velocity range with coincident image contours
(e.g., contours at a 3$\sigma$ level). We derive an integrated flux density of
11.24$\pm$0.09 Jy km s-1, 1.12$\pm$0.07 Jy km s-1, and 0.45$\pm$0.07 Jy km
s-1, by integrating the intensity over the line-emitting region in each
channel of the spectra shown in Figure 4, for F07599C, F07599S, and F07599NE,
respectively. From this analysis, the two subcomponents (F07599S and F07599NE)
are identified at 16.0$\sigma$ and 6.4$\sigma$ in the intensity maps of
residual data cube, respectively. These weak features are also clearly present
in the channel maps, where extended structures are found in the S and NE
direction (see Appendix B).

Figure 5: Left: optical images of IRAS F07599+6508. From left to right and top
to bottom, we show the SDSS $u$-, $g$-, $r$-, $i$-, $z$-band, and Pan-STARRS1
$y$-band images. Right: overlay of CO(1$-$0) contours (colored) of the three
components on the HST/WFPC2 F702W band image. The CO(1$-$0) contours are at
the same levels as in Figure 4. The optical diffraction spikes in the middle
of the object reflect the central luminous point source of AGN that is
significantly brighter than its host galaxy. Both the SDSS $gri$ and HST
optical images show clear signs of extended structure to the southwest and
southeast with respect to the main source.

Figure 5 shows the multi-wavelength images of IRAS F07599+6508 in SDSS $u$-,
$g$-, $r$-, $i$-, $z$-, and PanSTARRS1 $y$-bands, and an overlay of the
CO(1$-$0) emission of the three decomposed components on the HST/WFPC2 F702W
optical image. The optical morphologies shown in $gri$ and HST images are
found to be complex with tidal-like features, indicative of an asymmetric
stellar distribution, while the ultraviolet (UV; $u$-band at $\sim$ 0.35
$\mu$m) and near-IR ($z$\- and $y$-bands at $\sim$ 1 $\mu$m) emission show
more regular and compact distribution. A consistent morphology structure in
the UV band has been presented by Surace & Sanders (2000). The tidal features
seen in the optical is probably connected with an ongoing or past merging
event in IRAS F07599+6508. Moreover, the absence of tidal tails observed in
the near-IR bands is likely to suggest that the emission from optical tidal
tails may trace the massive, young stellar population formed at recent epochs,
since the near-IR emission is mostly sensitive to the longer-lived, low-mass
stars. If this is indeed the case, it will also explain the morphology
observed in UV band where the emission is largely attenuated in dust obscured
regions.

Interestingly, by comparing the distribution of CO molecular gas and stellar
emission (see the right panel of Figure 5), we find that the tidal feature in
the optical is coincident with the CO substructure F07599S. This provides
evidence that the CO feature observed in the S region may be a real structure.
Moreover, the CO emission also show an extension to the southeast (SE), the
same direction as the other tidal-tail-like feature observed in the optical
images, although the CO extended structure toward the SE appear in only a few
adjacent channels (e.g., channel maps at the velocity of about $-$178 km s-1
and $-$58 km s-1; see Appendix B). Meanwhile, a faint peak coincident with the
location of the CO peak in F07599NE is identified in the HST image, indicative
of a possible association.

####  3.2.3 CO kinematics

In Figure 6 we show the velocity and velocity dispersion maps of CO(1$-$0)
emission within [$-$330,+330] km s-1 velocity range for the whole galaxy IRAS
F07599+6508 and the three decomposed components, respectively. These maps were
created by applying a flux threshold of 3$\sigma$. The velocity map of the
whole galaxy IRAS F07599+6508 shows irregular structure, indicative of complex
velocity distribution. For the main component F07599C that is decomposed from
the whole system, a systematic velocity gradient seems to be present in the
north-south direction, from about $-$70 to +70 km s-1, although the gas at the
edge to the southeast and west appear to be disturbed. This implies that the
bulk of molecular gas is likely rotationally supported. However, the spatial
resolution of our data is not adequate to unambiguously interpret the observed
velocity gradient as due to disk rotation. An alternative possibility could be
that a velocity gradient and broad line profile may originate from an ongoing
merger, as the velocity gradient is observed roughly along the direction
connecting F07599C and F07599S. For F07599S, a velocity gradient of $\sim$180
km s-1(from about $-$70 to +110 km s-1) is also observed in the north-south
direction, while no well-defined structure is visible in the velocity map of
F07599NE, in which the velocity is positive around 60$-$90 km s-1 with respect
to the QSO redshift.

Figure 6: CO(1$-$0) velocity field (top) and velocity dispersion (bottom) maps
for the whole galaxy IRAS F07599+6508 (first column) and the three
subcomponents, F07599C (second column), F07599S (third column), and F07599NE
(fourth column). The CO velocity and velocity dispersion maps are calculated
using pixels deteccted at $\geqslant 3\sigma$ and are masked where the
integrated intensity is $<3\sigma$ in each component. The plus sign marks the
peak of 3 mm continuum emission, while the two crosses indicate the CO peak of
F07599S and F07599NE, respectively. The synthesized beam is shown in the
bottom-left corner of the first column panel.

If we assume a rotating disk geometry for the main component F07599C, we can
estimate the dynamical mass within the CO-emitting region as $M_{\rm
dyn}/M_{\odot}\approx 2.3\times 10^{5}v^{2}_{\rm cir}R$ (e.g., Wang et al.,
2010), where $R$ is the CO disk radius of 4.6 kpc (0.75 times the CO FWHM
major axis from visibility fitting; see Table 1) and $v_{\rm cir}$ is the
maximum circular velocity of the gas disk in km s-1. The $v_{\rm cir}$ is
estimated as $v_{\rm cir}=0.75\Delta v_{\rm FWHM}/{\rm sin}i$ (Ho, 2007),
where $i$ is the inclination angle between the polar axis of disk and the line
of sight and can be roughly calculated from the CO minor and major axis ratio
($a_{\rm min}/a_{\rm maj}$; see Table 1), i.e., $i={\rm cos^{-1}}(a_{\rm
min}/a_{\rm maj})=45^{\circ}$. The derived $M_{\rm dyn}$ is $2.1\times
10^{11}\ M_{\odot}$ with a $v_{\rm cir}$ of 450 km s-1 for the primary disk
F07599C. We stress that this estimate of dynamical mass can be significantly
affected by the large uncertainties of source size, inclination angle, and the
assumption of an inclined disk geometry. The total amount of molecular gas in
the host galaxy is $M_{\rm tot,mol}\sim 1.1\times 10^{10}\ M_{\odot}$ by
adopting a classic ULIRG-like CO-to-H2 conversion factor of $\alpha_{\rm
CO}=0.8\ M_{\odot}\ {\rm(K\ km\ s^{-1}\ pc^{-2})^{-1}}$. We obtain a molecular
gas fraction of $\sim 0.05$ for IRAS F07599+6508 by assuming the gas fraction
can be approximated as the ratio of gas to dynamical mass, slightly lower than
those (median value of 0.22$\pm$0.04) of IR QSOs measured by Tan et al.
(2019), and the local ULIRGs as well (Downes & Solomon, 1998). Here we note
the large uncertainties involved in the estimate of molecular gas mass, due to
the assumption of conversion factor $\alpha_{\rm CO}$, which is uncertain and
varies significantly from source to source (Bolatto et al., 2013). Recent
studies show evidence for higher $\alpha_{\rm CO}$ values in (U)LIRGs that are
affected by outflows if a substantial amount of molecular gas is in dense,
gravitationally bound states (e.g., Cicone et al., 2018, and references
therein).

The velocity dispersion map shows $\sigma_{\rm v}\sim$160 km s-1 at the center
of IRAS F07599+6508, which is significantly higher than the typical peak
dispersion of $\sim$80 km s-1 observed in local IR QSOs (Tan et al., 2019),
and $z\gtrsim$6 QSOs ($\sim$130 km s-1) as well (Neeleman et al., 2021). A
plausible explanation could be related to the beam smearing effect that
increases the observed dispersion, since the spatial angular resolution of our
NOEMA observations is about an order of magnitude lower than the ALMA
observations ($\sim 0.\arcsec 5$) shown by Tan et al. (2019). Assuming a
uniform intrinsic dispersion across the galaxy disk of IRAS F07599+6508, we
then estimate a dispersion of $\sigma_{\rm v}\sim 70-100$ km s-1, which is
measured in the outer region of F07599C where the beam smearing is considered
to be less severe (Davies et al., 2011). The velocity dispersion map of
F07599S shows an increase in the north direction with a peak dispersion of
$\sim$110 km s-1 which is in the overlap region between F07599C and F07599S,
while the gas in F07599NE appears less perturbed with dispersion of $\sim
30-60$ km s-1. The increased velocity dispersion observed in the overlap
region could be interpreted as a turbulence-dominated shock region induced by
a merger or a possible CO outflow component. Higher resolution observations
are needed to better constrain the intrinsic velocity dispersion.

###  3.3 Radio Continuum Emission

In star-forming galaxies, the radio emission is dominated by synchrotron
radiation from relativistic electrons accelerated in supernovae remnants of
massive stars. It is largely optically thin and unaffected by dust extinction,
and therefore an excellent probe of very recent star formation activity in
galaxies (Condon, 1992). For the radio-quite (RQ) AGNs, which represent the
majority of the AGN population, the origin of radio emission could be
connected with a wide range of possible mechanisms, namely, star formation,
AGN-driven wind, jets, and accretion-disk coronae (Panessa et al., 2019).

Figure 7 shows the VLA 1.49 GHz contour maps overlaid on the NOEMA CO(1-0)
data that is decomposed into three components. Along the direction of the
minor axis of synthesized beam of 1.49 GHz observations, two clumpy structures
outside the central bright compact radio source are potentially detected at
$\sim 3\sigma$ and $\sim 4\sigma$ level to the S and NE direction,
respectively. The radio peaks of both the central compact source and the two
clumps are well aligned with the CO peaks of the three decomposed components,
reinforcing the fact that the faint CO features observed in the S and NE
direction are very likely to be real.

Figure 7: VLA 1.49 GHz contour (black) maps of IRAS F07599+6508 overlaid on
the CO(1$-$0) contour (colored) maps of the three components. The CO(1$-$0)
contours are at the same levels as in Figure 4. The contours of 1.49 GHz
continuum emission start at 2$\sigma$ (1$\sigma$ noise level is 74.9 $\mu$Jy
beam-1) and increase by a factor of 1.5. The solid orange and black ellipses
in the bottom left and right corners are the beam size of CO and 1.49 GHz
maps, with angular resolutions of 5.″61$\times$3.″77 and 5.″33$\times$2.″08,
respectively.

A 2D Gaussian fit to the 1.49 GHz image returns a peak flux density of
32.33$\pm$0.10 mJy beam-1 and a consistent integrated flux density of
32.43$\pm$0.19 mJy, indicative of a relatively compact structure of radio
emission. This measurement agrees well with the 1.4 GHz flux density of 39.5
mJy derived from the NRAO VLA Sky Survey (Condon et al., 1998). IRAS
F07599+6508 is additionally detected at 4.86, 8.44, 14.94, and 22.46 GHz by
VLA. The maps and the flux measurements at these radio frequencies are shown
in Appendix C. At the redshift $z=0.1486$, observations at the radio
frequencies listed in Table 2 sample $\sim$1.7, 5.6, 9.7, 17.2, and 25.8 GHz
in the rest-frame, respectively. In the following section, we adopt the rest-
frame frequencies for clarity.

Figure 8 shows the SED of the continuum emission at radio and mm bands in IRAS
F07599+6508. At GHz frequencies, the flux density $S_{\nu}$ in terms of the
observed frequency $\nu$ is generally defined as $S_{\nu}\propto\nu^{\alpha}$,
with $\alpha$ the radio spectral index. The spectral indices measured from
1.7$-$5.6, 5.6$-$9.7, 9.7$-$17.2, and 17.2$-$25.8 GHz are
$\alpha_{5.6}^{1.7}=-1.10\pm 0.09$, $\alpha_{9.7}^{5.6}=-1.16\pm 0.22$,
$\alpha_{17.2}^{9.7}=-1.18\pm 0.26$, and $\alpha_{25.8}^{17.2}=-0.34\pm 0.78$,
respectively. It is clear that the spectrum turns over from steep to flat at
17.2$-$25.8 GHz, which is consistent with the expectation that the
contribution of thermal emission increases at higher frequencies. At low
frequencies, the spectral index $\alpha_{5.6}^{1.7}$ is slightly flatter than
those at higher frequencies, $\alpha_{9.7}^{5.6}$ and $\alpha_{17.2}^{9.7}$,
both of which show almost identical values. The flatter index observed at
lower frequencies has been interpreted as a result of optical depth effects
that the radio spectrum typically shows a sharp decline due to free-free
absorption (e.g., Condon, 1992; Clemens et al., 2010). In contrast, the flux
at larger frequency should be much less affected (e.g., Clemens et al., 2008).
Therefore, we assume a synchrotron (non-thermal) spectral index of
$\alpha^{\rm NT}=-1.18$, equivalent to the steepest part of the radio spectrum
measured between 9.7 and 17.2 GHz (see the dashed line in Figure 8). This
slope is slightly steeper than a normal synchrotron spectral index of
$\alpha=-0.8$ derived for a normal star-forming galaxy (Condon, 1992), but
consistent with the slopes measured for RQ AGN (e.g., Behar et al., 2015).

We extrapolate the contribution of non-thermal emission from
$\alpha_{17.2}^{9.7}$ at the NOEMA frequency of 90.5 GHz ($\sim$ 103.9 GHz in
the rest-frame) and subtract from the observed flux density. In addition,
subtracting off an estimate of the thermal dust contribution to the 103.9 GHz,
which is derived from a grey-body model fit to the $Hershel$/SPIRE data with
an emissivity spectral index $\beta$=1.6 (Clements et al., 2018), we obtain an
estimate of free-free emission component of 0.17 mJy. This implies a free-free
radio fraction at 103.9 GHz of $\sim$30%. A similar procedure applied to 25.8
GHz data leads to a free-free fraction of $\sim$30% estimated at this
frequency. This is consistent with the median thermal fraction of 20-30% at 22
GHz measured for local ULIRGs (Vega et al., 2008).

We derive an extrapolated 1.4 GHz flux density of 43.4 mJy from the observed
rest-frame 1.7 GHz flux density by assuming a radio power-law index of -1.18.
Using the flux densities measured at IRAS 60 and 100 $\mu$m, we calculate a
FIR-radio flux ratio444The FIR-radio flux ratio $q$ is defined as the ratio of
the FIR luminosity (40-120 $\mu$m) to the radio power (Condon, 1992):
$q\equiv$ log${}_{10}\left(\frac{L_{\rm FIR}}{3.75\times 10^{12}\
\rm{W}}\right)$-log${}_{10}\left(\frac{L_{\rm 1.4GHz}}{\rm{W\
Hz^{-1}}}\right)$ $q=1.67$. This implies a linear ratio about 4 times smaller
than the average value of $q=2.34\pm 0.01$ found for local star-forming
galaxies (Yun et al., 2001), possibly suggestive of an additional contribution
from compact AGN source that boosting the radio flux.

Figure 8: The radio-to-mm spectral energy distribution of IRAS F07599+6508
based on the data from VLA (circle) archive and NOEMA (square).

##  4 Discussion

The main results of our NOEMA observations are the detection of 3 mm continuum
emission and the multiple spatially decomposed CO components in the host
galaxy of IRAS F07599+6508. By combining with multi-wavelength data, we begin
the discussion by looking into the properties and the possible origin of these
components.

###  4.1 Nature of the 3mm Emission

The 3 mm continuum emission of IRAS F07599+6508 is resolved with a deconvolved
size of 1.″5 ($\sim 4.0$ kpc). Assuming a non-thermal spectral index of
$\alpha^{\rm NT}=-1.18$, we find that about 30% of the emission at 3 mm arises
from a thermal free-free component. This is significantly lower than the
thermal fraction of $\sim 70\%$ at 3 mm expected for normal star-forming
galaxies (Condon, 1992), indicating that the nonthermal synchrotron component
from AGN may contribute substantially to the 3 mm emission in IRAS
F07599+6508. This could be expected since it does not follow the FIR-radio
correlation characteristic of star-forming galaxies, with a ratio of $L_{\rm
FIR}/L_{\rm 1.4GHz}$ being about four times lower. A similar situation has
been observed in Mrk 231, the nearest IR QSO with a compact radio core(e.g.,
Lonsdale et al., 2003), where a rather low thermal fraction is observed at 33
GHz (Barcos-Muñoz et al., 2017). As mentioned above, an additional nonthermal
component from AGN could contribute to the emission at high radio frequencies.
An alternative explanation for the observed low thermal fraction is connected
with the possibility that, the production of free-free emission might be
suppressed due to the absorption of ionizing stellar photons by dust that
concentrates in the dense starburst regions.

Taking the 2-10 keV X-ray luminosity from Brightman & Nandra (2011), we find a
radio ($\sim 3$ mm) luminosity to X-ray luminosity ratio of $L_{\rm mm}/L_{\rm
X}\sim 3\times 10^{-2}$ for IRAS F07599+6508. This value is about two orders
of magnitude lower than those of radio-loud (RL) AGNs, which typically have
powerful relativistic jets that shine in the radio band, but close to the
ratios ($10^{-2}-10^{-4}$) measured for RQ AGNs (Behar et al., 2015; Panessa
et al., 2019), suggestive of a different origin of radio emission compared
with RL QSOs. We note that the ratio of $L_{\rm mm}/L_{\rm X}\sim 3\times
10^{-2}$ derived here could be an overestimated value, given that the AGN
embedded in the galaxy IRAS F07599+6508 is known to be heavily obscured and
thus a portion of X-ray photons may be absorbed (see Luo et al., 2014; La
Caria et al., 2019). A connection between the radio emission and the nuclear
accretion disk has been proposed to explain the correlation of the radio
luminosity $L_{\rm mm}$ and the X-ray luminosity $L_{\rm X}$ observed in RQ
AGNs (e.g., Laor & Behar, 2008). More evidences for an X-ray corona (i.e., hot
plasma with a temperature of $\sim 10^{9}$ K) origin of radio emission in RQ
AGNs have been presented by the detections of millimeter excess (e.g., Behar
et al., 2015; Inoue & Doi, 2018). These are in line with the theoretical
prediction for a radio synchrotron component from X-ray corona in the vicinity
of black hole accretion disks contributing to the emission in the millimeter
band (see Inoue & Doi, 2014; Raginski & Laor, 2016). Nevertheless, a
conclusive evidence for the coronal radio emission would be the detection of
the radio/millimeter and X-ray variability correlation similar to that
observed in stellar coronae (Panessa et al., 2019).

The thermal free-free emission arises directly from Hii regions ionized by
massive stars and its intensity is proportional to the ionizing photon rate
from young star-forming regions, and thereby provide an excellent diagnostic
for the current SFR of galaxies (Condon, 1992). Taking the free-free emission
flux measured at rest-frame 103.9 GHz and the free-free SFR calibration given
in Murphy et al. (2012), we calculate an SFR of 77 $M_{\odot}$ yr-1 for IRAS
F07599+6508 by assuming a thermal spectral index of $-$0.1 and an electron
temperature of $T_{\rm e}=10^{4}$ K. Compared to the IR-derived SFR of 280
$M_{\odot}$ yr-1 that is estimated based on the heated dust component (Hao et
al., 2008), a lower SFR measured by the free-free emission could be expected
since free-free emission is sensitive to massive stars with ages $\lesssim$10
Myr and consequently probe the very recent star formation activity. In
addition, the free-free-derived SFR would be underestimated if a significant
fraction of ionizing photons are absorbed by dust. Alternatively, it could be
that a portion of IR luminosity arise from the dust heating by nuclear
luminous AGN, which would lead to an overestimate of SFR. Indeed, we obtain an
SFR of 150 $M_{\odot}$ yr-1 using the calibration given by Murphy et al.
(2012) by correcting for the AGN contribution to the bolometric luminosity
($\alpha_{\rm AGN}=0.75$; Veilleux et al., 2009), which is close to the free-
free SFR estimate. All of these including the uncertainties (e.g., variations
in the actual $T_{\rm e}$) on the estimation of free-free-derived SFR may
explain the discrepancy in SFR estimates.

###  4.2 Origin of the off-center CO components

By exploiting the CO line cube of IRAS F07599+6508, we have identified two
off-center CO subcomponents, F07599S and F07599NE, at a significance level of
16.0$\sigma$ and 6.4$\sigma$ in the intensity map, respectively. Figure 4
shows the CO(1$-$0) spectra and maps extracted for these two subcomponents and
the main component F07599C as well. F07599S and F07599NE are located to the
south and northeast of the disk component with a projected distance of
$\sim$11.4 and 19.1 kpc and with a fraction of 9% and 4% of CO flux
contributed to the total CO flux, respectively. These two CO subcomponents are
not detected in 3mm continuum emission. For the main component F07599C, the
systematic CO velocity gradient observed may be a signature of a rotating disk
structure. If confirmed, the source size of $\sim 6.1$ kpc measured would
suggest that the CO disk is very extended, which is in sharp contrast with
local ULIRGs where the majority of molecular gas is centrally concentrated
within a radius of $\sim$1 kpc (Downes & Solomon, 1998; Wilson et al., 2008).

The optical images of IRAS F07599+6508 show evidence for a tidal-tail-like
structure (see Figure 5), suggestive of a recent merging event in IRAS
F07599+6508. It is interesting to note that the tidal feature shown in the
optical image is spatially coincident with the CO feature detected in the S
region. In addition, a faint peak coincident with the location of the CO peak
in F07599NE is identified in the HST optical image, likely associated with the
origin of the CO plume feature observed in the NE direction. It can be seen
from the channel maps in Appendix B that the CO emission is extended to the S
and NE direction, suggesting that the CO emission identified in the
extranuclear region of F07599S and F07599NE are very likely to be real
components. Furthermore, the two potentially detected clumpy structures
revealed in the 1.49 GHz radio image are found to align well with the peaks of
F07599S and F07599NE (see Figure 7), as well as the faint substructures
detected at 4.86 and 8.44 GHz (see Appendix C), supporting for the reality of
such CO feature.

In Figure 4 we have shown that the CO(1$-$0) spectrum extracted for F07599S is
best-fitted with a two-component (narrow and broad) Gaussian. The detection of
narrow (FWHM = 100$\pm$10 km s-1) feature may indicate the presence of
extended dynamical “quiescent” CO gas. For the broad component, the Gaussian
fit gives an FWHM of 500$\pm$80 km s-1, slightly larger than the line FWHM
(430$\pm$10 km s-1) of the primary component F07599C. Together with the shift
in velocity ($\sim$100 km s-1) and the tidal-tail-like structure observed in
the optical image, it is likely that F07599S represents a separate component
from main galaxy, and has a merger origin. One of the possibilities is that
F07599S may trace the emission from a small obscured companion to the host
galaxy of IRAS F07599+6508. The non-detection of UV emission in F07599S show
evidence for the large attenuation by dust (see Figure 5). More evidence for a
merger origin is given by the increased velocity dispersion observed in the
overlap region between F07599S and the main component, which could be
interpreted as a turbulence-dominated shock region induced by an interaction.

For F07599NE, the best-fit of the CO(1$-$0) line is a single Gaussian with
both FWHM and peak velocity consistent within the errors with the main
component. The significant lower velocity dispersion observed in this
component imply that F07599NE may have a different origin from F07599S.
Similar gas substructure has been observed in high-$z$ starbursts in [C ii]
emission but with a much smaller physical separation ($\sim 2$ kpc; Tadaki et
al., 2020), where the formation of subcomponent is expected to be associated
with a gravitationally unstable disk or shocks. By contrast, the large
physical separation (projected distance of $\sim 19$ kpc) observed between
F07599NE and the galaxy nucleus indicates that the CO plume feature observed
in the NE direction could be a separate component, likely associated with the
tidal debris reminiscent of a previous merger. Similar results have been seen
in the IR QSO IRAS F23411+0228, where an off-center clumpy structure is
detected in both CO line and 3 mm continuum maps and is found to align with a
faint peak identified in the HST image, which was interpreted as either an
ongoing merger or a clumpy star-forming region at large radii (Tan et al.,
2019). Cicone et al. (2015) have shown a very extended (up to $r\sim$ 30 kpc)
[C ii] cold gas component observed in a $z>6$ QSO host galaxy and they
interpret that the extended [C ii] emission may arise from a different ISM
phase than dense photon-dominated regions.

We note that there is a source accompanying IRAS F07599+6508 in multi-
wavelength data set (from GALEX near-UV to WISE 3.4 $\mu$m bands; see Appendix
D), with a distance of 17.″5 (projected distance of $\sim$ 45 kpc) to the
northwest. This source has a cross-identification as WISEA J080428.28+650001.7
(hereafter J0804+6500) in the NASA/IPAC555http://ned.ipac.caltech.edu
Extragalactic Database. There is no optical spectroscopic data available for
this object. In Appendix D, we show the NOEMA spectrum extracted from the
position of J0804+6500 by fitting the visibilites with a point source. We do
not detect any significant emission in either continuum or line from this
source. With the data available, it is not possible to assess whether
J0804+6500 is physically related with the IRAS F07599+6508 system, or is only
a projected companion source. The detection of line emission from this source
would help us measure the redshift and further determine whether it is an
interacting companion of IRAS F07599+6508.

Together with the CO clumpy structures of local IR QSOs revealed with ALMA
high resolution observations (Tan et al., 2019), the multiple components
revealed in the CO maps show evidence for the complex morphology that may
associate with the intense star formation and AGN activities. However,
additional deep CO observations with higher angular resolution are needed for
confirming the faint features observed in F07599S and F07599NE and testing the
corresponding hypothesis.

Figure 9: Left: continuum-subtracted spectra of the CO(1$-$0) emission line
(top) and the residuals (bottom) after removing a Gaussian fit (red line in
the top panel) to the CO spectrum of IRAS F07599+6508. In the top panel, the
NOEMA CO spectrum (black) is obtained by summing up the flux of the three
decomposed components in each channel (see the left panel in Figure 4), while
the IRAM 30m CO spectrum (blue) is also plotted for a comparison. In the
bottom panel, the dashed lines denote the 1$\sigma$ level. The spectra are
binned to a velocity resolution of $\sim$60 km s-1. Right: contours from
CO(1$-$0) dirty maps integrated over the velocity range of [$-$1340,$-$2000]
km s-1. Contours starts at $\pm 2\sigma$ and increase in steps of 1$\sigma$
(81 $\mu$Jy beam-1) with positive (negative) contours shown as solid (dotted)
lines. The plus sign and cross indicate the CO peak positions of F07599C and
F07599S, respectively. The synthesized beam is shown in the bottom-left
corner.

###  4.3 Comparison of molecular gas with Mrk 231

As noted previously in Section 1, IRAS F07599+6508 and Mrk 231 share many
common or similar properties in multi-wavelength, i.e., X-ray, UV, optical,
and IR, and have been proposed as local templates of candidate ULIRG-to-
optical QSO transition objects. High-resolution CO observations of Mrk 231
show that the molecular gas is concentrated in the central $\sim$1 kpc
rotating disk (Downes & Solomon, 1998; Cicone et al., 2012; Feruglio et al.,
2015), though O2 emission might be present over $\sim$10 kpc (Wang et al.,
2020). Compared with Mrk 231, the molecular gas of IRAS F07599+6508 is also
likely distributed in a disk with systematic velocity gradient in the main
component, but much more extended with a size of $\sim$6.1 kpc. In addition,
the CO line width of IRAS F07599+6508 is $\sim$430 km s-1 FWHM, much broader
than that of 190 km s-1 measured for the core component of Mrk 231 (Downes &
Solomon, 1998). The CO(2$-$1) to CO(1$-$0) line luminosity ratios derived for
IRAS F07599+6508 and Mrk 231 are both $\sim 0.8$ (Xia et al., 2012; Cicone et
al., 2012), indicative of modest sub-thermal molecular gas in their hosts.

The detection of massive, highly energetic, and kpc-scale molecular outflows
in Mrk 231 show direct observational evidence for AGN feedback onto the host
galaxy (e.g., Feruglio et al., 2010, 2015; Aalto et al., 2012; Cicone et al.,
2012), consistent with the AGN-galaxy co-evolutionary models which propose
that the outflow feedback is an important process in quenching the star
formation activity in massive galaxies (e.g., Di Matteo et al., 2005; Hopkins
et al., 2006). The CO(1$-$0) wing-to-peak ratio in Mrk 231 is $\sim$2.5%, in
agreement with the ratio of 1-8% typically measured for local ULIRGs (Feruglio
et al., 2010; Cicone et al., 2012, 2014).

The presence of molecular outflow in IRAS F07599+6508 has been confirmed by
the detection of OH 119 $\mu$m P-Cygni profile with velocities $<$ 1000 km s-1
(Veilleux et al., 2013). A fast wind has also been seen in this system in the
ionized gas as traced by the broad blueshift optical spectra (H$\alpha$ and
H$\beta$) and in the neutral atomic phase with broad absorption observed in Na
i D line (e.g., Zheng et al., 2002; Rupke et al., 2017).

In the left-top panel of Figure 9 we show the continuum-subtracted NOEMA
CO(1-0) spectrum of the whole galaxy of IRAS F07599+6508, which is obtained by
summing up the flux of the three decomposed components in each channel,
compared with the CO(1-0) line observed with the IRAM 30m (Xia et al., 2012).
The CO peak flux measured by NOEMA is found to agree well with the IRAM 30m
data, suggesting that little emission is filtered out by our interferometric
observations. The residual spectrum after removing a Gaussian fit to the core
component of the CO spectrum is shown in the left-bottom panel of Figure 9. A
noise level (1$\sigma$) of 0.68 mJy in each $\sim$ 60 km s-1 channel was
measured, corresponding to 2.3% of the peak flux density. Scaling to the
signal strengths from Mrk 231 would lead to an estimated flux density of the
CO broad emission of $\sim$ 0.75 mJy for IRAS F07599+6508, comparable to the
noise level of our NOEMA data.

It is interesting to note that several consecutive channels in the blue-
shifted side of the CO spectrum (see Figure 9, left-bottom) show flux density
higher than 1$\sigma$. We averaged the visibilities of the high-velocity
component in the velocity range [$-$1340, $-$2000] km s-1 and fitted a point-
source model. The peak of high-velocity CO emission is located between the
F07599C and F07599S with a projected distance of 4.″3 from the galaxy nucleus
(see Figure 9, right). The flux derived from uv_fit is 0.228$\pm$0.054 Jy km
s-1, corresponding to a detection significance of 4.2$\sigma$. However, this
does not satisfy the criterion defined by Cicone et al. (2014), i.e., CO
emission with velocities $>500$ km s-1 is detected with a significance of at
least 5$\sigma$, for claiming the detection of a molecular outflow. The OH
outflow counterpart in CO emission is not significantly detected in our NOEMA
data, although there is a hint of some flux at very high blueshifted
velocities in the overlap region between F07599C and F07599S (see Figure 9),
where the extracted CO spectrum shows a clear decomposition into broad plus
narrow components (see the CO spectrum of F07599S in Figure 4). Nevertheless,
it is difficult to rule out the existence of a weak CO outflow with the
current data.

Given the much extended CO disk providing fuel for active large-scale star
formation in IRAS F07599+6508, perhaps a much weaker outflow feedback in
quenching star formation is expected here as compared to that of Mrk 231,
where bursts of star formation only remain in the central compact core. If
that is true, IRAS F07599+6508 should be in an earlier phase than that of Mrk
231 in the sense of ULIRGs evolution in transition into the classical QSOs. We
stress that additional NOEMA CO observations are required to investigate this
source in detail, particularly the presence of a possible CO outflow.

##  5 Summary

In this work, we present the deep NOEMA interferometric observations of
CO(1$-$0) line and 3 mm continuum emission in the ultraluminous IR QSO IRAS
F07599+6508 at $z=0.1486$. This system shares many physical properties with
Mrk 231, the nearest and most luminous IR QSO, in multi-wavelength. Both the
CO line and 3mm continuum emission are resolved in our observations, and the
CO gas is found to be spatially decomposed into three components. Combining
these observations with results from multi-wavelength data (i.e., radio,
optical, and infrared), we study the properties of the 3 mm continuum emission
and the CO components spatially decomposed, and the main results are
summarized as follows.

1\. The 3 mm continuum emission is resolved with a deconvolved size of 1.″5
($\sim 4.0$ kpc). Assuming a non-thermal spectral index of $\alpha^{\rm
NT}=-1.18$, equal to the slope measured between 9.7 and 17.2 GHz representing
the steepest part of the radio spectrum, we find that $\sim 30$% of the
emission at 3 mm arises from a free-free component, significantly lower than
the fraction ($\sim 70\%$) expected for a normal star-forming galaxy. This
suggests that the nonthermal synchrotron component from AGN may contribute
substantially to the 3 mm emission. The relatively low FIR-to-radio flux ratio
of $q=1.67$ derived also show evidence for an additional contribution from the
central AGN that boosting the radio flux. We estimate an SFR of 77 $M_{\odot}\
{\rm yr^{-1}}$ based on the thermal free-free emission component at 3 mm,
which is close to the IR estimate of $150\ M_{\odot}\ {\rm yr^{-1}}$ corrected
for the AGN contribution to the bolometric luminosity.

2\. Two off-center CO subcomponents, F07599S and F07599NE, located to the
south and northeast of the main component F07599C with a projected distance of
$\sim$11.4 and 19.1 kpc, respectively, were identified by model fits to the
$uv$-data. More evidence for the reality of F07599S and F07599NE are given by
the spatial coincidence of the tidal-tail-like feature and the faint peak in
the optical images, and the potentially detected clumpy structures in the 1.49
GHz radio image with these two components. The fraction of F07599S and
F07599NE contributed to the total CO flux are 9% and 4%, respectively.
Analysis of the CO kinematics together with the multi-wavelengths data for the
decomposed components suggests that the gas in the dominant component F07599C
is likely rotationally supported with a systematic velocity gradient observed,
while F07599S may represent a separate component associated with a merger and
F07599NE probably have a different origin.

3\. Careful analysis of our current CO(1$-$0) data indicates that, the
significance of the high-velocity emission feature seen in the blue-shifted
side of the spectrum is not sufficient for claiming the detection of a CO
outflow in IRAS F07599+6508. It should be noted that, however, the spatially
broad CO core distribution with a size larger than 6 kpc is in marked contrast
with the compact $\sim$ kpc CO concentration in Mrk 231, where strong
molecular outflows were observed over wide velocity ranges.

Additional deeper observations with higher angular resolution are required to
better constrain the CO morphology and kinematics on kpc and sub-kpc scales,
and confirm the presence of high-velocity gas and further determine whether
the features seen here are potential signs of outflowing gas that is escaping
from the host galaxy.

We thank the anonymous referee for the constructive comments that helped
improve the paper. We acknowledge the support from IRAM staff during the
observations and data reduction. This work is based on observations carried
out under project number S18CA with the IRAM NOEMA Interferometer. IRAM is
supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain). This work was
supported by National Key Basic Research and Development (R&D) Program of
China (Grant No. 2017YFA0402704), NSFC Grant No. 11803090, 11861131007,
12033004, and 11420101002, and Chinese Academy of Sciences Key Research
Program of Frontier Sciences (Grant No. QYZDJ-SSW-SLH008). X.Y.X and C.N.H.
acknowledge support from the NSFC (Grant No. 11733002). The work by K.K is
supported by the JSPS KAKENHI Grant Number JP17H06130 and the NAOJ ALMA
Scientific Research Grant Number 2017-06B.

##  Appendix A Decomposition of CO emission

To illustrate the spatial decomposition of CO emission in IRAS F07599+6508, we
here describe the process of model fitting procedure. We fit the model to the
$uv$-data using combinations of functions by allowing all parameters to be
free. Combinations of models in a successive progression of complexity started
with a point source and then circular Gaussian, and elliptical Gaussian. The
simultaneous fit was rerun each time as necessary to achieve a uniform
residual (i.e., no significant peak or negative regions $\gtrsim 3\sigma$)
after subtraction of the models. By verifying the evolution of the
uncertainties on the fitted parameters, we found that the best fit was
obtained with a combination of an elliptical Gaussian (F07599C) and two point
sources (F07599S and F07599NE) functions with minimum uncertainties on the
fitted parameters. The residuals after subtracting an elliptical Gaussian plus
two point source components are shown in the right panel of Figure 3. It
should be noted that there are a few peaks with SNR$\sim 3$ (e.g., the peaks
to the northwest of F07599S and F07599NE, respectively) in the residual map,
indicative of the likely presence of extended structure. These weak features
are however difficult to be tested with the current data. We should stress
that higher spatial resolution observations are needed to confirm these
findings.

##  Appendix B Channel maps of CO(1$-$0) emission

The spatial and velocity distribution of CO(1$-$0) emission in IRAS
F07599+6508 can be further explored by inspecting the channel maps in Figure
10, which shows the image of CO emission within an interval of $\sim$120 km
s-1 encompassing the velocity range from $-$2090 km s-1 to 2093 km s-1. The
channel maps were produced from the continuum-subtracted data cube that was
cleaned down to the 1$\sigma$ noise level. The extended structures to the
south and northeast of galaxy center observed in each channel map, spatially
coincident with the subcomponent of F07599S and F07599NE (marked as red
crosses in Figure 10), show evidence for the reality of such features.

##  Appendix C Radio maps and flux measurements

Table 2 summarizes the radio continuum properties of IRAS F07599+6508 obtained
with VLA. Except the 22.46 GHz data that was obtained in a single CD
configuration, observations at the rest of frequencies were made with multiple
arrays.

To investigate the distribution of radio emission and make a comparison with
CO emission, we show the map obtained with higher angular resolution for each
frequency in Figure 11. For 14.94 and 22.46 GHz, the emission are detected
with a compact morphology, while for the emission at 4.86 and 8.44 GHz where
the observations were deeper, the radio source exhibits more complex
morphology with clumpy structures outside the main body of the galaxy to the
east and west, respectively. The peak flux density measured in different
configuration were found to be in good agreement for each frequency (i.e., to
within $\sim 10\%$; see Table 2). The highest resolution image was obtained at
17.2 GHz with beam size of $\sim$0.″18 ($\sim$ 0.47 kpc), of which the peak
flux density measured is found to be comparable with that measured from CD
configuration (beam size $\sim$4.″58). The consistent peak flux densities
measured between different beams for the different bands indicate that the
bulk of radio emission probably originates in the inner compact region of the
AGN.

Figure 10: Channel maps showing the CO(1$-$0) emission of IRAS F07599+6508 in
the [$-$2000, +1700] km s-1 range with channel wide of $\sim$120 km s-1. The
red crosses mark the positions of CO subcomponents (F07599S and F07599NE) and
the black plus sign indicates the peak of 3 mm continuum emission. Contour
levels start at $\pm 1.5\sigma$ and increase by factors of 1.5 with positive
(negative) contours shown as solid (blue dotted) lines. The velocities in km
s-1 are marked in the top left corner of each channel map and the synthesized
beam is shown in the bottom-left corner.

Figure 11: VLA multi-band radio continuum imaging of IRAS F07599+6508 at 4.86, 8.44, 14.94, and 22.46 GHz from left to right. Except the observations at 22.46 GHz which was conducted in a single array of CD configuration, the rest of observations were made with arrays in both AB configuration (A configuration for 14.94 GHz) and CD configuration. Here we only show the maps with higher angular resolution. All the contours start at 3$\sigma$ and increase by a factor of 1.5. The synthesized beam and beam size of each radio observations are shown in the bottom of each panel. The white plus sign indicates the peak of 3 mm continuum emission. Table 2: VLA radio continuum properties of IRAS F07599+6508  Obs. freq. | RAaaSource position derived from a 2D Gaussian fit (using imfit task in CASA) to the high-resolution images. | DECaaSource position derived from a 2D Gaussian fit (using imfit task in CASA) to the high-resolution images. | Conf. | Synt.beam | SbbPeak flux density measured from imfit using a single 2D Gaussian component.peak
---|---|---|---|---|---
(GHz) | (J2000) | (J2000) |  | (arcsec) | (mJy beam-1)
1.49 GHz | 08:04:30.53 | 64:59:53.3 | AB | 5.33$\times$2.08 | 32.33$\pm$0.10
|  |  | CD | 65.2$\times$28.3 | 36.72$\pm$0.16
4.86 GHz | 08:04:30.51 | 64:59:52.6 | AB | 1.84$\times$0.94 | 9.59$\pm$0.05
|  |  | CD | 20.8$\times$8.8 | 9.15$\pm$0.16
8.44 GHz | 08:04:30.51 | 64:59:52.6 | AB | 1.10$\times$0.47 | 4.67$\pm$0.07
|  |  | CD | 11.0$\times$5.1 | 5.23$\pm$0.09
14.94 GHz | 08:04:30.471 | 64:59:52.76 | A | 0.18$\times$0.13 | 2.38$\pm$0.06
|  |  | CD | 4.58$\times$2.56 | 2.69$\pm$0.08
22.46 GHz | 08:04:30.46 | 64:59:52.8 | CD | 2.92$\times$1.93 | 2.20$\pm$0.10

Figure 12: Image cutouts (60″$\times$60″on a side) of IRAS F07599+6508 and a
companion source WISEA J080428.28+650001.7 to the northwest. From left to
right, we show the GALEX near-UV, SDSS $u$-, $r$-, Pan-STARRS1 $y$-, 2MASS
Ks-, and WISE 3.4$\mu$m-band images.

Figure 13: NOEMA spectrum (without correction for PBA) of WISEA
J080428.28+650001.7 covered in the USB (top) and the LSB (bottom). The spectra
are extracted by fitting a point-source model to the $uv$-data and are binned
to a channel width of $\sim 120-140$ km s-1.

##  Appendix D Multi-wavelength data of IRAS F07599+6508 and a companion
source to the northwest

In this appendix we show the multi-wavelength (from GALEX near-UV to WISE
3.4$\mu$m band) data of IRAS F07599+6508 and a companion source WISEA
J080428.28+650001.7 to the northwest (Figure 12). J0804+6500 has a distance of
17.″5 ($\sim$ 45 kpc) from IRAS F07599+6508, with a coordinates of
08:04:28.140, +65:00:02.25.

Figure 13 shows the NOEMA spectra (without correction for primary beam
attenuation (PBA)) extracted from the position of J0804+6500 by fitting a
point-source model to the $uv$-data. The tuning frequency covered in the upper
and lower sidebands corresponds to a redshift range of 0.128 to 0.223 and
0.238 to 0.462, respectively, if associated with a CO(1$-$0) transition. The
1$\sigma$ rms noise level of the spectrum in the USB and LSB shown in Figure
13 is 0.20 mJy and 0.21 mJy, respectively. We do not detect any significant
emission from J0804+6500 in our observations.

## References

  * Aalto et al. (2012) Aalto, S., Garcia-Burillo, S., Muller, S., et al. 2012, A&A, 537, A44, doi: 10.1051/0004-6361/201117919
  * Barcos-Muñoz et al. (2017) Barcos-Muñoz, L., Leroy, A. K., Evans, A. S., et al. 2017, ApJ, 843, 117, doi: 10.3847/1538-4357/aa789a
  * Barnes & Hernquist (1992) Barnes, J. E., & Hernquist, L. 1992, ARA&A, 30, 705, doi: 10.1146/annurev.aa.30.090192.003421
  * Behar et al. (2015) Behar, E., Baldi, R. D., Laor, A., et al. 2015, MNRAS, 451, 517, doi: 10.1093/mnras/stv988
  * Bolatto et al. (2013) Bolatto, A. D., Wolfire, M., & Leroy, A. K. 2013, ARA&A, 51, 207, doi: 10.1146/annurev-astro-082812-140944
  * Boroson & Meyers (1992) Boroson, T. A., & Meyers, K. A. 1992, ApJ, 397, 442, doi: 10.1086/171800
  * Brightman & Nandra (2011) Brightman, M., & Nandra, K. 2011, MNRAS, 414, 3084, doi: 10.1111/j.1365-2966.2011.18612.x
  * Bryant & Scoville (1996) Bryant, P. M., & Scoville, N. Z. 1996, ApJ, 457, 678, doi: 10.1086/176763
  * Canalizo & Stockton (2001) Canalizo, G., & Stockton, A. 2001, ApJ, 555, 719, doi: 10.1086/321520
  * Cao et al. (2008) Cao, C., Xia, X. Y., Wu, H., et al. 2008, MNRAS, 390, 336, doi: 10.1111/j.1365-2966.2008.13747.x
  * Cicone et al. (2012) Cicone, C., Feruglio, C., Maiolino, R., et al. 2012, A&A, 543, A99, doi: 10.1051/0004-6361/201218793
  * Cicone et al. (2020) Cicone, C., Maiolino, R., Aalto, S., Muller, S., & Feruglio, C. 2020, A&A, 633, A163, doi: 10.1051/0004-6361/201936800
  * Cicone et al. (2014) Cicone, C., Maiolino, R., Sturm, E., et al. 2014, A&A, 562, A21, doi: 10.1051/0004-6361/201322464
  * Cicone et al. (2015) Cicone, C., Maiolino, R., Gallerani, S., et al. 2015, A&A, 574, A14, doi: 10.1051/0004-6361/201424980
  * Cicone et al. (2018) Cicone, C., Severgnini, P., Papadopoulos, P. P., et al. 2018, ApJ, 863, 143, doi: 10.3847/1538-4357/aad32a
  * Clemens et al. (2010) Clemens, M. S., Scaife, A., Vega, O., & Bressan, A. 2010, MNRAS, 405, 887, doi: 10.1111/j.1365-2966.2010.16534.x
  * Clemens et al. (2008) Clemens, M. S., Vega, O., Bressan, A., et al. 2008, A&A, 477, 95, doi: 10.1051/0004-6361:20077224
  * Clements et al. (2018) Clements, D. L., Pearson, C., Farrah, D., et al. 2018, MNRAS, 475, 2097, doi: 10.1093/mnras/stx3227
  * Condon (1992) Condon, J. J. 1992, ARA&A, 30, 575, doi: 10.1146/annurev.aa.30.090192.003043
  * Condon et al. (1998) Condon, J. J., Cotton, W. D., Greisen, E. W., et al. 1998, AJ, 115, 1693, doi: 10.1086/300337
  * Daddi et al. (2010) Daddi, E., Bournaud, F., Walter, F., et al. 2010, ApJ, 713, 686, doi: 10.1088/0004-637X/713/1/686
  * Davies et al. (2011) Davies, R., Förster Schreiber, N. M., Cresci, G., et al. 2011, ApJ, 741, 69, doi: 10.1088/0004-637X/741/2/69
  * Di Matteo et al. (2005) Di Matteo, T., Springel, V., & Hernquist, L. 2005, Nature, 433, 604, doi: 10.1038/nature03335
  * Downes & Solomon (1998) Downes, D., & Solomon, P. M. 1998, ApJ, 507, 615, doi: 10.1086/306339
  * Evans et al. (2001) Evans, A. S., Frayer, D. T., Surace, J. A., & Sanders, D. B. 2001, AJ, 121, 1893, doi: 10.1086/319972
  * Evans et al. (2006) Evans, A. S., Solomon, P. M., Tacconi, L. J., Vavilkin, T., & Downes, D. 2006, AJ, 132, 2398, doi: 10.1086/508416
  * Fabian (2012) Fabian, A. C. 2012, ARA&A, 50, 455, doi: 10.1146/annurev-astro-081811-125521
  * Feruglio et al. (2010) Feruglio, C., Maiolino, R., Piconcelli, E., et al. 2010, A&A, 518, L155, doi: 10.1051/0004-6361/201015164
  * Feruglio et al. (2015) Feruglio, C., Fiore, F., Carniani, S., et al. 2015, A&A, 583, A99, doi: 10.1051/0004-6361/201526020
  * Fischer et al. (2010) Fischer, J., Sturm, E., González-Alfonso, E., et al. 2010, A&A, 518, L41, doi: 10.1051/0004-6361/201014676
  * Greve et al. (2005) Greve, T. R., Bertoldi, F., Smail, I., et al. 2005, MNRAS, 359, 1165, doi: 10.1111/j.1365-2966.2005.08979.x
  * Guilloteau & Lucas (2000) Guilloteau, S., & Lucas, R. 2000, in Astronomical Society of the Pacific Conference Series, Vol. 217, Imaging at Radio through Submillimeter Wavelengths, ed. J. G. Mangum & S. J. E. Radford, 299
  * Hao et al. (2005) Hao, C. N., Xia, X. Y., Mao, S., Wu, H., & Deng, Z. G. 2005, ApJ, 625, 78, doi: 10.1086/429716
  * Hao et al. (2008) Hao, C.-N., Xia, X.-Y., Shu-DeMao, Deng, Z.-G., & Wu, H. 2008, Chinese J. Astron. Astrophys., 8, 12, doi: 10.1088/1009-9271/8/1/02
  * Ho (2007) Ho, L. C. 2007, ApJ, 668, 94, doi: 10.1086/521270
  * Hopkins et al. (2006) Hopkins, P. F., Hernquist, L., Cox, T. J., et al. 2006, ApJS, 163, 1, doi: 10.1086/499298
  * Hopkins et al. (2008) Hopkins, P. F., Hernquist, L., Cox, T. J., & Kereš, D. 2008, ApJS, 175, 356, doi: 10.1086/524362
  * Imanishi et al. (2019) Imanishi, M., Nakanishi, K., & Izumi, T. 2019, ApJS, 241, 19, doi: 10.3847/1538-4365/ab05b9
  * Inoue & Doi (2014) Inoue, Y., & Doi, A. 2014, PASJ, 66, L8, doi: 10.1093/pasj/psu079
  * Inoue & Doi (2018) —. 2018, ApJ, 869, 114, doi: 10.3847/1538-4357/aaeb95
  * Jarvis et al. (2020) Jarvis, M. E., Harrison, C. M., Mainieri, V., et al. 2020, MNRAS, 498, 1560, doi: 10.1093/mnras/staa2196
  * Kim & Sanders (1998) Kim, D. C., & Sanders, D. B. 1998, ApJS, 119, 41, doi: 10.1086/313148
  * King & Pounds (2015) King, A., & Pounds, K. 2015, ARA&A, 53, 115, doi: 10.1146/annurev-astro-082214-122316
  * Kormendy & Ho (2013) Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511, doi: 10.1146/annurev-astro-082708-101811
  * Kwan et al. (1995) Kwan, J., Cheng, F.-Z., Fang, L.-Z., Zheng, W., & Ge, J. 1995, ApJ, 440, 628, doi: 10.1086/175302
  * La Caria et al. (2019) La Caria, M. M., Vignali, C., Lanzuisi, G., Gruppioni, C., & Pozzi, F. 2019, MNRAS, 487, 1662, doi: 10.1093/mnras/stz1381
  * Laor & Behar (2008) Laor, A., & Behar, E. 2008, MNRAS, 390, 847, doi: 10.1111/j.1365-2966.2008.13806.x
  * Lawrence et al. (1999) Lawrence, A., Rowan-Robinson, M., Ellis, R. S., et al. 1999, MNRAS, 308, 897, doi: 10.1046/j.1365-8711.1999.02593.x
  * Le Floc’h et al. (2005) Le Floc’h, E., Papovich, C., Dole, H., et al. 2005, ApJ, 632, 169, doi: 10.1086/432789
  * Lonsdale et al. (2003) Lonsdale, C. J., Lonsdale, C. J., Smith, H. E., & Diamond, P. J. 2003, ApJ, 592, 804, doi: 10.1086/375778
  * Luo et al. (2014) Luo, B., Brandt, W. N., Alexander, D. M., et al. 2014, ApJ, 794, 70, doi: 10.1088/0004-637X/794/1/70
  * McMullin et al. (2007) McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, in Astronomical Society of the Pacific Conference Series, Vol. 376, Astronomical Data Analysis Software and Systems XVI, ed. R. A. Shaw, F. Hill, & D. J. Bell, 127
  * Murphy et al. (2012) Murphy, E. J., Bremseth, J., Mason, B. S., et al. 2012, ApJ, 761, 97, doi: 10.1088/0004-637X/761/2/97
  * Neeleman et al. (2021) Neeleman, M., Novak, M., Venemans, B. P., et al. 2021, arXiv e-prints, arXiv:2102.05679.  https://arxiv.org/abs/2102.05679
  * Panessa et al. (2019) Panessa, F., Baldi, R. D., Laor, A., et al. 2019, Nature Astronomy, 3, 387, doi: 10.1038/s41550-019-0765-4
  * Raginski & Laor (2016) Raginski, I., & Laor, A. 2016, MNRAS, 459, 2082, doi: 10.1093/mnras/stw772
  * Rupke et al. (2017) Rupke, D. S. N., Gültekin, K., & Veilleux, S. 2017, ApJ, 850, 40, doi: 10.3847/1538-4357/aa94d1
  * Rupke & Veilleux (2011) Rupke, D. S. N., & Veilleux, S. 2011, ApJ, 729, L27, doi: 10.1088/2041-8205/729/2/L27
  * Sanders & Mirabel (1996) Sanders, D. B., & Mirabel, I. F. 1996, ARA&A, 34, 749, doi: 10.1146/annurev.astro.34.1.749
  * Sanders et al. (1988a) Sanders, D. B., Soifer, B. T., Elias, J. H., et al. 1988a, ApJ, 325, 74, doi: 10.1086/165983
  * Sanders et al. (1988b) Sanders, D. B., Soifer, B. T., Elias, J. H., Neugebauer, G., & Matthews, K. 1988b, ApJ, 328, L35, doi: 10.1086/185155
  * Scoville et al. (2003) Scoville, N. Z., Frayer, D. T., Schinnerer, E., & Christopher, M. 2003, ApJ, 585, L105, doi: 10.1086/374544
  * Shangguan et al. (2020) Shangguan, J., Ho, L. C., Bauer, F. E., Wang, R., & Treister, E. 2020, ApJS, 247, 15, doi: 10.3847/1538-4365/ab5db2
  * Silk & Rees (1998) Silk, J., & Rees, M. J. 1998, A&A, 331, L1.  https://arxiv.org/abs/astro-ph/9801013
  * Soifer et al. (1987) Soifer, B. T., Sanders, D. B., Madore, B. F., et al. 1987, ApJ, 320, 238, doi: 10.1086/165536
  * Solomon et al. (1997) Solomon, P. M., Downes, D., Radford, S. J. E., & Barrett, J. W. 1997, ApJ, 478, 144, doi: 10.1086/303765
  * Spoon et al. (2013) Spoon, H. W. W., Farrah, D., Lebouteiller, V., et al. 2013, ApJ, 775, 127, doi: 10.1088/0004-637X/775/2/127
  * Surace & Sanders (2000) Surace, J. A., & Sanders, D. B. 2000, AJ, 120, 604, doi: 10.1086/301491
  * Tadaki et al. (2020) Tadaki, K.-i., Iono, D., Yun, M. S., et al. 2020, ApJ, 889, 141, doi: 10.3847/1538-4357/ab64f4
  * Tan et al. (2019) Tan, Q.-H., Gao, Y., Kohno, K., et al. 2019, ApJ, 887, 24, doi: 10.3847/1538-4357/ab50be
  * Taylor et al. (1999) Taylor, G. B., Silver, C. S., Ulvestad, J. S., & Carilli, C. L. 1999, ApJ, 519, 185, doi: 10.1086/307362
  * Ueda et al. (2014) Ueda, J., Iono, D., Yun, M. S., et al. 2014, ApJS, 214, 1, doi: 10.1088/0067-0049/214/1/1
  * van Moorsel et al. (1996) van Moorsel, G., Kemball, A., & Greisen, E. 1996, in Astronomical Society of the Pacific Conference Series, Vol. 101, Astronomical Data Analysis Software and Systems V, ed. G. H. Jacoby & J. Barnes, 37
  * Vega et al. (2008) Vega, O., Clemens, M. S., Bressan, A., et al. 2008, A&A, 484, 631, doi: 10.1051/0004-6361:20078883
  * Veilleux et al. (1999) Veilleux, S., Kim, D. C., & Sanders, D. B. 1999, ApJ, 522, 113, doi: 10.1086/307634
  * Veilleux et al. (2009) Veilleux, S., Rupke, D. S. N., Kim, D. C., et al. 2009, ApJS, 182, 628, doi: 10.1088/0067-0049/182/2/628
  * Veilleux et al. (2013) Veilleux, S., Meléndez, M., Sturm, E., et al. 2013, ApJ, 776, 27, doi: 10.1088/0004-637X/776/1/27
  * Véron-Cetty et al. (2006) Véron-Cetty, M. P., Joly, M., Véron, P., et al. 2006, A&A, 451, 851, doi: 10.1051/0004-6361:20054237
  * Wang et al. (2020) Wang, J., Li, D., Goldsmith, P. F., et al. 2020, ApJ, 889, 129, doi: 10.3847/1538-4357/ab612d
  * Wang et al. (2010) Wang, R., Carilli, C. L., Neri, R., et al. 2010, ApJ, 714, 699, doi: 10.1088/0004-637X/714/1/699
  * Wilson et al. (2008) Wilson, C. D., Petitpas, G. R., Iono, D., et al. 2008, ApJS, 178, 189, doi: 10.1086/590910
  * Xia et al. (2001) Xia, X.-Y., Boller, T., Deng, Z.-G., & Börner, G. 2001, Chinese J. Astron. Astrophys., 1, 221, doi: 10.1088/1009-9271/1/3/221
  * Xia et al. (2012) Xia, X. Y., Gao, Y., Hao, C.-N., et al. 2012, ApJ, 750, 92, doi: 10.1088/0004-637X/750/2/92
  * Xu et al. (2014) Xu, C. K., Cao, C., Lu, N., et al. 2014, ApJ, 787, 48, doi: 10.1088/0004-637X/787/1/48
  * Xu et al. (2015) —. 2015, ApJ, 799, 11, doi: 10.1088/0004-637X/799/1/11
  * Yun et al. (2001) Yun, M. S., Reddy, N. A., & Condon, J. J. 2001, ApJ, 554, 803, doi: 10.1086/323145
  * Zheng et al. (2002) Zheng, X. Z., Xia, X. Y., Mao, S., Wu, H., & Deng, Z. G. 2002, AJ, 124, 18, doi: 10.1086/340964