# ALMA CO Observations of the Mixed-Morphology Supernova Remnant W49B:
Efficient Production of Recombining Plasma and Hadronic Gamma-rays via
Shock–Cloud Interactions

H. Sano National Astronomical Observatory of Japan, Mitaka, Tokyo 181-8588,
Japan<EMAIL_ADDRESS>S. Yoshiike Department of Physics, Nagoya
University, Furo-cho, Chikusa-ku, Nagoya 464-8601, Japan  Y. Yamane Department
of Physics, Nagoya University, Furo-cho, Chikusa-ku, Nagoya 464-8601, Japan
K. Hayashi Institute of Space and Astronautical Science (ISAS), Japan
Aerospace Exploration Agency (JAXA), 3-1-1, Yoshinodai, Chuo-ku, Sagamihara,
Kanagawa 252-5210, Japan  R. Enokiya Department of Physics, Institute of
Science and Technology, Keio University, 3-14-1 Hiyoshi, Kohoku-ku, Yokohama,
Kanagawa 223-8522, Japan  K. Tokuda National Astronomical Observatory of
Japan, Mitaka, Tokyo 181-8588, Japan<EMAIL_ADDRESS>Department of
Physical Science, Graduate School of Science, Osaka Prefecture University, 1-1
Gakuen-cho, Naka-ku, Sakai 599-8531, Japan  K. Tachihara Department of
Physics, Nagoya University, Furo-cho, Chikusa-ku, Nagoya 464-8601, Japan  G.
Rowell School of Physical Sciences, The University of Adelaide, North Terrace,
Adelaide, SA 5005, Australia  M. D. Filipović Western Sydney University,
Locked Bag 1797, Penrith South DC, NSW 1797, Australia  Y. Fukui Department of
Physics, Nagoya University, Furo-cho, Chikusa-ku, Nagoya 464-8601, Japan

###### Abstract

We carried out new CO($J$ = 2–1) observations toward the mixed-morphology
supernova remnant (SNR) W49B with the Atacama Large Millimeter/submillimeter
Array (ALMA). We found that CO clouds at $\sim$10 km s-1 show a good spatial
correspondence with synchrotron radio continuum as well as an X-ray deformed
shell. The bulk mass of molecular clouds accounts for the western part of the
shell, not for the eastern shell where near-infrared H2 emission is detected.
The molecular clouds at $\sim$10 km s-1 show higher kinetic temperature of
$\sim$20–60 K, suggesting that modest shock-heating occurred. The expanding
motion of the clouds with $\Delta V\sim$6 km s-1 was formed by strong winds
from the progenitor system. We argue that the barrel-like structure of Fe rich
ejecta was possibly formed not only by an asymmetric explosion, but also by
interactions with dense molecular clouds. We also found a negative correlation
between the CO intensity and the electron temperature of recombining plasma,
implying that the origin of the high-temperature recombining plasma in W49B
can be understood as the thermal conduction model. The total energy of
accelerated cosmic-ray protons $W_{\mathrm{p}}$ is estimated to be
$\sim$$2\times 10^{49}$ erg by adopting an averaged gas density of
$\sim$$650\pm 200$ cm-3. The SNR age–$W_{\mathrm{p}}$ diagram indicates that
W49B shows one of the highest in-situ values of $W_{\mathrm{p}}$ in the gamma-
ray bright SNRs.

Supernova remnants (1667); Interstellar medium (847); Cosmic ray sources
(328); Gamma-ray sources (633); X-ray sources (1822)

††software: CASA (v 5.6.0.: McMullin et al., 2007), CIAO (v 4.12: Fruscione et
al., 2006), CALDB (v 4.9.1 Graessle et al., 2007).††facilities: ALMA, Chandra,
NuSTAR, Very Large Array (VLA), Nobeyama 45-m Telescope, James Clerk Maxwell
Telescope (JCMT).

##  1 INTRODUCTION

Mixed-morphology supernova remnants (SNRs) are characterized a center-filled
thermal X-ray morphology with a synchrotron radio shell, accounting for more
than 25% of Galactic SNRs (Rho & Petre, 1998). Most of the mixed-morphology
SNRs are interacting with the dense interstellar medium (ISM) evidenced by
radio-line emission such as CO, H I, and/or 1720 MHz OH masers (e.g., Seta et
al., 1998; Yusef-Zadeh et al., 2003; Kuriki et al., 2018). In addition, some
of them are associated with GeV/TeV gamma-ray sources, which likely arise from
interactions between accelerated cosmic-ray (CR) protons and dense clouds in
the vicinity of the SNRs (e.g., Aharonian et al., 2008; Bamba et al., 2016).
Moreover, shock-propagation into the clumpy ISM and/or dense circumstellar
matter (CSM) have a potential to explain their mixed-morphology and thermal
X-ray radiation (e.g., Shimizu et al., 2012; Slavin et al., 2017; Zhang et
al., 2019). Therefore, the shock-interacting ISM plays an important role in
understanding their morphology, plasma conditions, and cosmic-ray acceleration
(see also reviews by Vink, 2012; Yamaguchi, 2020; Sano & Fukui, 2021). To
unveil the physical processes and high-energy phenomena in the mixed-
morphology SNRs, detailed comparative studies among the radio-line emission,
X-rays, and gamma-rays are needed.

W49B (also known as G43.3$-$0.2) is a well-studied Galactic mixed-morphology
SNR with the bright radio-continuum shell and thermal-dominated center-filled
X-rays as shown in Figure 1. The SNR is thought to be lying on the far-side of
the Galaxy from us (e.g., Lockhart & Goss, 1978; Brogan & Troland, 2001). The
small apparent diameter of $\sim$$3^{\prime}$–$5^{\prime}$ is consistent with
the larger distance of $\sim$7.5–11.3 kpc (Zhu et al., 2014; Ranasinghe &
Leahy, 2018; Lee et al., 2020) and its young age (5–6 kyr, Zhou & Vink, 2018).
W49B is also thought to be an efficient accelerator of CR protons owing to its
bright GeV/TeV gamma-rays with the pion-decay bump (H. E. S. S. Collaboration
et al., 2018). The total energy of CR protons was derived to be
$\sim$1049–1051 erg, assuming the targeted gas density of 10–1000 cm-3 (Abdo
et al., 2010).

Figure 1: Intensity map of Chandra broad-band X-rays ($E$: 0.5–7.0 keV, e.g.,
Lopez et al., 2013b) superposed on the VLA radio continuum contours at 20 cm
obtained from MAGPIS (Helfand et al., 2006). The contour levels are 5.0, 7.6,
15.4, 28.4, 46.6, and 70.0 mJy beam-1. The region enclosed by dashed rectangle
corresponds to the observed region with ALMA ACA.

The X-ray radiation of W49B is characterized by three properties: the most
luminous in Fe K-shell line emission (Yamaguchi et al., 2014), non-thermal
Bremsstrahlung (Tanaka et al., 2018), and recombining (overionized) plasma
where the ionization temperature goes even higher than the electron
temperature (Ozawa et al., 2009; Miceli et al., 2010; Lopez et al., 2013a;
Zhou & Vink, 2018; Yamaguchi et al., 2018; Sun & Chen, 2020; Holland-Ashford
et al., 2020; Siegel et al., 2020). The elongated structure of Fe-rich ejecta
is believed to be related to a bipolar/jet-driven Type Ib/Ic explosion and/or
interactions between the shock and a surrounding interstellar cloud (Keohane
et al., 2007; Lopez et al., 2013b; Bear & Soker, 2017). On the other hand,
recent X-ray studies on metal abundances favor Type Ia models (Zhou & Vink,
2018; Siegel et al., 2020). Additionally, Sun & Chen (2020) conclude that the
SN type is unclear, with neither core-collapse or Ia models perfectly
reproducing their best-fit abundances. The origin of recombining
plasma—thermal conduction with cold-dense clouds and/or adiabatic cooling—is
still being debated (e.g., Yamaguchi et al., 2018; Sun & Chen, 2020; Holland-
Ashford et al., 2020). If we detect decreasing the electron temperature toward
the shocked clouds, we can confirm the thermal conduction scenario as a
formation mechanism of the recombination plasma (e.g., Matsumura et al.,
2017b; Okon et al., 2018, 2020).

Although W49B is thought to be interacting with interstellar clouds, it is a
perplexing question which clouds are physically associated. Keohane et al.
(2007) discovered 2.12 $\mu$m shocked H2 emission toward the eastern and
southwestern shells by near-infrared photometric observations. The total mass
of shocked H2 is estimated to 14–550 $M_{\sun}$. Subsequent radio observations
using CO line emission revealed three molecular clouds at velocities of
$\sim$10, $\sim$40, and $\sim$60 km s-1, which are possibly associated with
the SNR (e.g., Zhu et al., 2014; Kilpatrick et al., 2016; Lee et al., 2020).
Zhu et al. (2014) argued that the molecular cloud at $\sim$40 km s-1 is
interacting with W49B because of its wind-bubble like morphology. On the other
hand, Kilpatrick et al. (2016) found a line-broadening of 12CO profile in a
molecular cloud at $\sim$10 km s-1 located toward the western shell of W49B,
and hence the authors claimed that the cloud at 10 km s-1 is interacting with
the SNR. The velocity is roughly consistent with the H I absorption studies
(Brogan & Troland, 2001; Ranasinghe & Leahy, 2018). Most recently, Lee et al.
(2020) performed near-infrared spectroscopy of shocked H2 emission toward four
strips on W49B. The authors found that a central velocity of shocked H2 is
$\sim$64 km s-1 and then concluded that the molecular cloud at $\sim$60 km s-1
located toward the center and the southwest shell of W49B is likely associated
with the SNR. In either case, detailed spatial and kinematic studies as well
as deriving cloud properties (e.g., mass, density, kinetic temperature) have
not been performed due to the modest sensitivity and angular resolution of CO
datasets up to $\sim$$20^{\prime\prime}$, corresponding to a spatial
resolution of $\sim$1 pc at the distance of 10 kpc.

In the present paper, we report on results of new millimeter wavelength
observations using CO($J$ = 2–1) line emission with the Atacama Compact Array
(ACA, also known as Morita Array) which is a part of the Atacama Large
Millimeter/submillimeter Array (ALMA). The unprecedented sensitivity and high-
angular resolution of $\sim$$7^{\prime\prime}$ ($\sim$0.3 pc at the distance
of 10 kpc) of the ALMA CO data enable us to identify the interacting molecular
cloud and its physical relation to the high-energy phenomena in W49B. Section
2 describes the observational datasets and reductions. Section 3 comprises
five subsections: Sections 3.1–3.3 present overview distributions of X-ray,
radio continuum, and CO; Section 3.4–3.5 show physical conditions of molecular
clouds. Discussion and conclusions are given in Sections 4 and 5,
respectively.

##  2 OBSERVATIONS AND DATA REDUCTIONS

###  2.1 CO

Observations of 12CO($J$ = 2–1) and 13CO($J$ = 2–1) line emission were
conducted using ALMA ACA Band 6 (211–275 GHz) as a Cycle 6 project (proposal
no. 2018.1.01780.S). We used the mosaic observation mode with 10–12 antennas
of 7-m array and four antennas of 12-m total power (TP) array. The observed
areas were $5\farcm 1\times 2\farcm 7$ rectangular regions centered at
($\alpha_{\mathrm{J2000}}$, $\delta_{\mathrm{J2000}}$) $=$
($19^{\mathrm{h}}11^{\mathrm{m}}09\fs 00$, $+9\arcdeg 06\arcmin 24\farcs 8$)
and ($19^{\mathrm{h}}11^{\mathrm{m}}07\fs 44$, $+9\arcdeg 05\arcmin 56\farcs
2$). The actual observed area is shown in Figure 1. The combined baseline
length of 7-m array data is from 8.85 to 48.95 m, corresponding to $u$–$v$
distances from 6.8 to 37.6 $k\lambda$ at 230.538 GHz. Two quasars,
J1924$-$2914 and J1751$+$0939, were observed as bandpass and flux calibrators.
We also observed four quasars, J1907$+$0127, J1922$+$1530, J1938$+$0448, and
J1851$+$0035, as phase calibrators. We performed data reduction using the
Common Astronomy Software Application (CASA, McMullin et al., 2007) package
version 5.5.0. We utilized “tclean” task with multi-scale deconvolver and
natural weighting. The emission mask was also selected using the auto-
multithresh procedure (Kepley et al., 2020). We combined the cleaned 7-m array
data and the calibrated TP array data using “feather” task to recover the
missing flux and diffuse emission. The beam size of feathered data is $8\farcs
23\times 4\farcs 77$ with a position angle of $-75.31\arcdeg$ for the 12CO($J$
= 2–1) data, and $8\farcs 28\times 5\farcs 04$ with a position angle of
$-79.37\arcdeg$ for the 13CO($J$ = 2–1) data. The typical noise fluctuations
are $\sim$0.065 K for the 12CO($J$ = 2–1) data and $\sim$0.055 K for the
13CO($J$ = 2–1) data at the velocity resolution of 0.4 km s-1.

We also used archival datasets of 12CO($J$ = 1–0) and 12CO($J$ = 3–2) line
emission for estimating physical properties of molecular clouds. The 12CO($J$
= 1–0) data are from the FOREST Unbiased Galactic Plane Imaging survey with
the Nobeyama 45 m telescope (FUGIN, Umemoto et al., 2017), and the 12CO($J$ =
3–2) data are from the CO High-Resolution Survey (COHRS, Dempsey et al., 2013)
obtained with the James Clerk Maxwell Telescope (JCMT). The angular resolution
is $\sim$$20^{\prime\prime}$ for the 12CO($J$ = 1–0) data and
$\sim$16$\farcs$6 for the 12CO($J$ = 3–2) data. The velocity resolutions of
12CO($J$ = 1–0) and 12CO($J$ = 3–2) data are 1.3 and 1.0 km s-1, respectively.
To improve the signal-to-noise ratio of the 12CO($J$ = 3–2) data, we combined
four spatial pixels and a rebind pixel size is to be $12^{\prime\prime}$. The
typical noise fluctuations are $\sim$1.4 K for the 12CO($J$ = 1–0) data and
$\sim$0.18 K for the 12CO($J$ = 3–2) data at each velocity resolution.

###  2.2 Radio Continuum

The radio continuum data at 20 cm wavelength are from the Multi-Array Galactic
Plane Imaging Survey (MAGPIS, Helfand et al., 2006) obtained with the Very
Large Array (VLA) and the Effelsberg 100-m telescope. The angular resolution
is $\sim$$6^{\prime\prime}$, which is compatible for the ALMA ACA resolution.
The typical noise fluctuations are $\sim$1–2 mJy.

###  2.3 X-rays

We utilized archival X-ray data obtained by Chandra (observation IDs are 117,
13440, and 13441), which have been published in numerous papers (e.g.,
Kawasaki et al., 2005; Lopez et al., 2009a, b, 2011, 2013a, 2013b; Yang et
al., 2009; Koo et al., 2016; Zhou & Vink, 2018). The X-ray datasets were taken
with the Advanced CCD Imaging Spectrometer S-array (ACIS-S3). We used Chandra
Interactive Analysis of Observations (CIAO, Fruscione et al., 2006) software
version 4.12 with CALDB 4.9.1 (Graessle et al., 2007) for data reduction and
imaging. After reprocessing for all datasets using the “chandra_repro”
procedure, we created an energy-filtered, exposure-corrected image using the
“fluximage” procedure in the energy bands of 0.5–7.0 keV (broad-band, see
Figure 1), 0.5–1.2 keV (soft-band), 1.2–2.0 keV (medium band), 2.0–7.0 keV
(hard-band), and 4.2–5.5 keV (continuum band). Because the soft- and medium-
band images are heavily affected by interstellar absorption (e.g., Zhou &
Vink, 2018), in this paper we focus on the X-ray images at energies greater
than 2.0 keV. We also created an exposure-corrected, continuum-subtracted
image of Fe He$\alpha$ line emission (6.4–6.9 keV) following the method
presented by Lopez et al. (2013b). The typical angular resolution of Chandra
images is $\sim$0$\farcs$5.

Figure 2: Integrated intensity maps of ALMA ACA 12CO($J$ = 2–1) (upper panels)
and 13CO($J$ = 2–1) (lower panels) for (a, b) the 10 km s-1 cloud, (c, d) the
40 km s-1 cloud, and (e, f) the 60 km s-1 cloud. The integration velocity
range is from 1 to 15 km s-1 for the 10 km s-1 cloud; from 38 to 47 km s-1 for
the 40 km s-1 cloud; and is from 57 to 67 km s-1 for the 60 km s-1 cloud.
Superposed contours are the same as shown in Figure 1.

##  3 RESULTS

###  3.1 Overview of X-ray, Radio Continuum, and CO Distributions

Figure 1 shows the Chandra broad-band X-ray image of W49B superposed on the
VLA radio continuum contours at 20 cm wavelength. As presented in previous
studies, a barrel-shaped radio-continuum shell with several co-axis filaments
and the center-filled X-rays are seen (e.g., Keohane et al., 2007; Lopez et
al., 2009a, 2011). The X-ray elongated feature brighter than $\sim 5\times
10^{-7}$ photon cm-2 s-1 is roughly consistent with the spatial distribution
of Fe He$\alpha$ emission (see Figure 3 in Lopez et al. 2013b). We note that
overall distributions of X-rays and radio continuum are quite different
between the northeastern and southwestern halves: the shell boundary of the
northeastern half roughly coincides with each other, whereas the southwestern
shell of X-rays is significantly deformed compared to that of radio continuum.
In particular, the X-rays are dim around positions at
($\alpha_{\mathrm{J2000}}$, $\delta_{\mathrm{J2000}}$) $\sim$
($19^{\mathrm{h}}11^{\mathrm{m}}00\fs 0$, $+09\arcdeg 06\arcmin
00^{\prime\prime}$), ($19^{\mathrm{h}}11^{\mathrm{m}}03\fs 0$, $+09\arcdeg
05\arcmin 00^{\prime\prime}$), and ($19^{\mathrm{h}}11^{\mathrm{m}}08\fs 0$,
$+09\arcdeg 06\arcmin 06^{\prime\prime}$): the first two correspond to the
bright peaks of radio-continuum and the other is partially surrounded by radio
filaments. This trend is also seen in the 4.2–5.5 keV band image which is
mostly free from the interstellar absorption as well as line-emission.

Figure 3: Velocity channel maps of ALMA ACA 12CO($J$ = 2–1) (upper panels) and
13CO($J$ = 2–1) (lower panels) for each cloud. Each panel shows CO integrated
intensity distribution integrated over the velocity range from 1 to 15 km s-1
every 2.8 km s-1 for the 10 km s-1 cloud; from 38 to 48 km s-1 every 2 km s-1
for the 40 km s-1 cloud; and 57 to 67 km s-1 every 2 km s-1 for the 60 km s-1
cloud. Superposed contours are the same as shown in Figure 1. Figure 3:
Continued.

Figure 2 shows integrated intensity maps of 12CO($J$ = 2–1) and 13CO($J$ =
2–1) for three velocity ranges of 1–15 km s-1 (hereafter “10 km s-1 cloud”),
38–47 km s-1 (hereafter “40 km s-1 cloud”), and 57–67 km s-1 (hereafter “60 km
s-1 cloud”) as previously mentioned in several papers (e.g., Zhu et al., 2014;
Kilpatrick et al., 2016; Lee et al., 2020). The kinematic distance of
molecular cloud is $\sim$11 kpc for the 10 km s-1 cloud; $\sim$9 kpc for the
40 km s-1 cloud 111Although Zhu et al. (2014) suggested the distance of 40 km
s-1 cloud to be $\sim$10 kpc using a Galactic rotation curve model with
$R_{0}=8.5$ kpc and $\Theta_{0}=220$ km s-1 (Kerr & Lynden-Bell, 1986), we
adopt its distance to be $\sim$9 kpc using the latest Galactic parameters of
$R_{0}=7.92$ kpc and $\Theta_{0}=227$ km s-1 (VERA Collaboration et al.,
2020). Here $R_{0}$ is the distance from the Sun to the Galactic center and
$\Theta_{0}$ is the rotation velocity of the local standard of rest. We use
the latter values throughout the paper.; and $\sim$7 kpc for the 60 km s-1
cloud (Sofue et al., 2020). Note that there are no other CO clouds within the
velocity range from $-15.0$ to 92.6 km s-1, and hence we focus on the three
molecular clouds in the present paper.

In the 10 km s-1 cloud (Figures 2a and 2b), there is an intensity gradient
increasing from southeast to northwest. The radio continuum shows fairly good
spatial correspondence with molecular clouds in the northern inward protrusion
and along the north, northwest, and southwest rims. On the other hand, both
the 12CO and 13CO emission lines are faint in the southeastern shell where
shocked H2 emission is strongly detected (e.g., Keohane et al., 2007). Dense
clouds traced by 13CO emission are located not only outside the shell
boundary, but also inside the radio continuum shell.

In the 40 km s-1 cloud (Figures 2c and 2d), the 12CO emission has a relatively
uniform distribution rather than that of the 10 km s-1 cloud, but the 13CO
emission shows an intensity gradient increasing from northwest to southeast.
We note that the radio brightest shell in west shows lack of both the 12CO and
13CO emission lines. The bright 12CO emission and 13CO clumps are located both
toward the southwest of the SNR as well as inside the southeastern shell. The
southwestern CO clumps seem to be along the sharp edge of radio shell, whereas
the CO clumps inside the SNR show no significant spatial correlations with the
radio shell morphology. At this spatial coverage, we could not find the
bubble-like CO structure toward W49B as mentioned by Zhu et al. (2014).

In the 60 km s-1 cloud (Figures 2e and 2f), there are dense clouds across the
SNR from northeast to southwest with two bright CO peaks at
($\alpha_{\mathrm{J2000}}$, $\delta_{\mathrm{J2000}}$) $\sim$
($19^{\mathrm{h}}11^{\mathrm{m}}17\fs 0$, $+09\arcdeg 07\arcmin
30^{\prime\prime}$) and ($19^{\mathrm{h}}10^{\mathrm{m}}57\fs 5$, $+09\arcdeg
05\arcmin 07^{\prime\prime}$). The former contains two H II regions cataloged
by Urquhart et al. (2009), whereas the latter does not have any cataloged
objects. The diffuse CO emission inside the SNR appears to be spatially anti-
correlated with the radio continuum contours.

To derive masses of three molecular clouds, we used the following equations:

$\displaystyle M=m_{\mathrm{H}}\mu\Omega D^{2}\sum_{i}N_{i}(\mathrm{H}_{2}),$
(1) $\displaystyle N(\mathrm{H}_{2})=X\cdot W(\mathrm{CO}),$ (2)

where $m_{\mathrm{H}}$ is the mass of hydrogen, $\mu=2.8$ is the mean
molecular weight, $\Omega$ is the solid angle of each pixel, $D$ is the
distance to W49B, $N_{i}(\mathrm{H}_{2})$ is the column density of molecular
hydrogen for each pixel, $X$ is CO-to-H2 conversion factor of $2\times
10^{20}$ cm-2 (K km s${}^{-1})^{-1}$ (Bertsch et al., 1993), and $W$(CO) is
the velocity integrated intensity of 12CO($J$ = 1–0) emission line obtained
from the FUGIN data (Umemoto et al., 2017). We estimated the mass of molecular
cloud inside the radio shell to be $\sim$$4.1\times 10^{4}$ $M_{\sun}$ for the
60 km s-1 cloud and $\sim$$2.7\times 10^{4}$ $M_{\sun}$ for the other two
clouds, where we adopted the shell radius of 2$\farcm$5 centered at
($\alpha_{\mathrm{J2000}}$, $\delta_{\mathrm{J2000}}$) $=$
($19^{\mathrm{h}}11^{\mathrm{m}}07\fs 34$, $+09\arcdeg 06\arcmin 01\farcs 1$).
These values are roughly consistent with previously derived cloud masses using
the 13CO($J$ = 1–0) emission and the 13CO-to-H2 conversion factor (H. E. S. S.
Collaboration et al., 2018).

###  3.2 Detailed Spatial Comparison with the Radio Continuum Shell

Figure 3 shows velocity channel maps for each molecular cloud superposed on
the radio continuum contours. We find that dense 13CO clumps at a velocity
rang from 3.8 to 9.4 km s-1 are nicely along not only with the northern and
southern shell, but also with radio filaments inside the SNR at positions of
($\alpha_{\mathrm{J2000}}$, $\delta_{\mathrm{J2000}}$) $\sim$
($19^{\mathrm{h}}11^{\mathrm{m}}06\fs 4$, $+09\arcdeg 07\arcmin
03^{\prime\prime}$) and ($19^{\mathrm{h}}11^{\mathrm{m}}05\fs 5$, $+09\arcdeg
05\arcmin 56^{\prime\prime}$). We also find that both the 12CO and 13CO clouds
at the velocity range of 9.4–12.2 km s-1 show global anti-correlation with the
radio shell. In addition, 12CO emission at the velocity range of 12.2–15.0 km
s-1 shows a good spatial correspondence with the western shell especially for
the sharp edge of the southwestern rim. This velocity range is consistent with
that inferred from Kilpatrick et al. (2016). The southwestern CO clumps at
42.0–44.0 km s-1 seem to be along the sharp edge of radio shell especially
prominent in the 13CO line emission. Moreover, the 13CO clouds at 38.0–40.0 km
s-1 shows the good spatial correspondence with the southeastern half of the
radio continuum shell. Furthermore, the 63–65 km s-1 12CO map shows a lack of
CO emission along the eastern shell and shows bright CO emission in the gaps
between the western and northern radio contours. Although the other CO clouds
also appear to be overlapped with the radio continuum shell and filaments,
their spatial correspondence is not clear.

Figure 4: Position-velocity diagram of ALMA ACA 12CO($J$ = 2–1) (upper panels)
and 13CO($J$ = 2–1) (lower panels) for (a, b) the 10 km s-1 cloud, (c, d) the
40 km s-1 cloud, and (e, f) the 60 km s-1 cloud. The integration range of
Right Ascension is from 283$\fdg$76 to 287$\fdg$79\. Dashed curves and circles
delineate expanding gas motion (see the text). Vertical dashed lines indicate
integration velocity ranges for each cloud.

###  3.3 Position–Velocity Diagrams

Figure 4 shows position–velocity diagrams for the three molecular clouds. The
velocity distributions of 12CO and 13CO emission are similar to each other
except for a velocity at $\sim$40 km s-1, suggesting that the 12CO emission at
$\sim$40 km s-1 is subject to self-absorption due to optically thick
component. In fact, the 12CO emission at $\sim$40 km s-1 shows dip-like
feature, whereas the 13CO spectrum at the same velocity has strong line
emission. We also find incomplete and complete cavity-like structures in the
10 and 40 km s-1 clouds, respectively (see dashed curves in Figures 4a–4d).
The incomplete cavity (or arc-like distribution) in the 10 km s-1 cloud lies
from 5 km s-1 to 11 km s-1 with a velocity dispersion of a few km s-1. On the
other hand, the complete cavity in the 40 km s-1 cloud is clearly seen
especially for 12CO emission, whose velocity range is from 41 km s-1 to 47 km
s-1. Note that spatial extents of these cavities are roughly consistent with
the diameter of the radio continuum shell. By contrast, there is no clear
evidence for such cavity-like structure in the position–velocity diagram of
the 60 km s-1 cloud (see Figures 4e and 4f).

###  3.4 Intensity Ratio Maps

Figure 5 shows intensity ratio maps of 12CO($J$ = 3–2) / 13CO($J$ = 2–1)
(hereafter $R_{\mathrm{12CO32/13CO21}}$) toward the three molecular clouds. A
higher value of $R_{\mathrm{12CO32/13CO21}}$ tends to be observed in diffuse
warm gas thermalized by supernova shocks and/or stellar radiation assuming
that the abundance ratio of 12CO/13CO is constant within a molecular cloud
complex (cf. Bieging & Peters, 2011; Dell’Ova et al., 2020). We find that
high-intensity ratios of $R_{\mathrm{12CO32/13CO21}}$ $\sim$10 are distributed
toward the shell of W49B in the 10 km s-1 cloud. The southwestern edge of the
shell (hereafter SW-edge) shows the highest value of
$R_{\mathrm{12CO32/13CO21}}$ $\sim$20, where the radio continuum shell is
strongly deformed. On the other hand, the 40 km s-1 cloud shows no significant
enhancement of $R_{\mathrm{12CO32/13CO21}}$ toward the SNR shell. The
southeast shell in the 60 km s-1 cloud also shows higher values of
$R_{\mathrm{12CO32/13CO21}}$, while the regions with high-intensity ratios
continuously extend beyond the radio-shell boundary and may not be related to
W49B.

Figure 5: Intensity ratio maps of 12CO($J$ = 3–2) / 13CO($J$ = 2–1) for (a)
the 10 km s-1 cloud, (b) the 40 km s-1 cloud, and (c) the 60 km s-1 cloud.
Each data was smoothed to match the beam size 16$\farcs$6\. The velocity range
of each cloud and superposed solid contours are the same as shown in Figure 1.
Superposed dashed contours indicate 13CO($J$ = 2–1) integrated intensities,
whose lowest contour levels are 3 K km s-1 and contour intervals are 1.0 K km
s-1 for the 10 km s-1 cloud; 0.5 K km s-1 for the 40 km s-1 cloud; and 1.5 K
km s-1 for the 60 km s-1 cloud. The gray areas represent that the 12CO($J$ =
3–2) and/or 13CO($J$ = 2–1) data show the lower significance than
$\sim$7$\sigma$. Yellow crosses (Reference SW/NW, SW/NE-shell, central-
filament) and a dashed circle (SW-edge) discussed in Section 3.5 are
indicated. Green squares and triangles also indicate the positions of H II
regions and IRAS point sources, respectively (Beichman et al., 1988; Urquhart
et al., 2009).

###  3.5 Physical conditions of the molecular clouds

To reveal physical conditions for each cloud in detail, we performed the Large
Velocity Gradient (LVG) analysis (e.g., Goldreich & Kwan, 1974; Scoville &
Solomon, 1974). The LVG analysis can calculate the radiative transfer of
molecular line emission, assuming a spherically isotropic cloud with uniform
photon escape probability, temperature, and radial velocity gradient of
$dv/dr$. Here, $dv$ is the half-width half-maximum of CO line profiles and
$dr$ is a cloud radius. We selected six individual CO peaks for each cloud,
which are significantly detected both the 12CO and 13CO emission lines. Four
of them appear to be along the radio continuum shell or filaments (hereafter
refer to as “shell clouds”), and the others were selected as reference which
are located outside of the shell (hereafter refer to as “reference clouds”).
CO spectra toward each position are shown in Figure 6. We adopt $dv/dr={2.5}$
km s-1 pc-1 for the 10 km s-1 SW-shell and the 40 km s-1 SE-shell; 1.5 km s-1
pc-1 for the 60 km s-1 SW-edge, NW-shell, and Reference-SW; and 1.0 km s-1
pc-1 for the others. We also utilized the abundance ratio of [12CO/H2] =
$5\times 10^{-5}$ (Blake et al., 1987) and the isotope abundance ratio of
[12CO/13CO] = 49 (Langer & Penzias, 1990).

Figure 6: CO intensity profiles toward individual CO peaks in (a–f) the 10 km
s-1 cloud, (g–l) the 40 km s-1 cloud, and (m–r) the 60 km s-1 cloud. The CO
spectra enclosed by vertical lines represent individual CO peaks that are
focused on the present study. Each CO data was smoothed to match the beam size
16$\farcs$6 and the velocity resolution of 1 km s-1. In SW-edge spectra of the
10 km s-1 cloud, we combined four pixels around the highest intensity ratio of
12CO($J$ = 3–2) / 13CO($J$ = 2–1) as shown in Figure 5a. Figure 7: LVG results
on the number density of molecular hydrogen, $n$(H2), and the kinetic
temperature, $T_{\mathrm{kin}}$ for each cloud as shown in Figure 6. The red
lines and blue dash-dotted lines indicate the intensity ratios of 12CO($J$ =
3–2) / 12CO($J$ = 2–1) and 12CO($J$ = 3–2) / 13CO($J$ = 2–1), respectively.
The shaded areas in red and blue represent $1\sigma$ error ranges of each
intensity ratio. Yellow crosses indicate the best-fit values of $n$(H2) and
$T_{\mathrm{kin}}$ for each cloud. The results are summarized in Table 1.

Figure 7 shows the LVG results on the number density of molecular hydrogen,
$n$(H2), and the kinematic temperature, $T_{\mathrm{kin}}$, toward the six
positions for each cloud. The best-fit values of $n$(H2) and
$T_{\mathrm{kin}}$ are summarized in Table 1. In the 10 km s-1 cloud, we find
that the shell clouds show $T_{\mathrm{kin}}$ $\sim$20–60 K, which are
significantly higher than that of the reference clouds
($T_{\mathrm{kin}}={15}$ K). By contrast, all shell clouds in both the 40 km
s-1 and 60 km s-1 components show $T_{\mathrm{kin}}$ $\sim$10 K which are
roughly consistent with their reference clouds except for Reference-SE in the
60 km s-1 ($T_{\mathrm{kin}}={21}$ K). We also note that there is no relation
between the number density of molecular hydrogen and kinetic temperature for
each cloud.

##  4 DISCUSSION

###  4.1 Molecular Clouds Associated with W49B

Previous studies proposed three candidates of molecular clouds which are
interacting with the SNR W49B: namely the 10 km s-1 cloud, 40 km s-1 cloud,
and the 60 km s-1 cloud (Zhu et al., 2014; Kilpatrick et al., 2016; Lee et
al., 2020). Their claim is mainly based on three elements: (1) a line-
broadening feature of CO emission at $\sim$10 km s-1, (2) a wind-bubble like
morphology of CO cloud at $\sim$40 km s-1, and (3) a central velocity of
shocked H2 line emission at $\sim$64 km s-1. In this section, we discuss which
cloud is the most likely to be associated with W49B in terms of spatial
distributions, presence of the expanding gas motion, and physical conditions
of CO clouds.

####  4.1.1 Spatial Distributions of CO Clouds

We first emphasize that the 10 km s-1 cloud shows a clear spatial
correspondence with the radio continuum shell and filaments (see Figure 3 and
Section 3.2). In particular, a majority of CO clouds at $\sim$10 km s-1 are
located along the outer boundary of the radio continuum shell: 13CO clouds in
the northern shell at 6.6–9.4 km s-1, and arc-like 12CO clouds in the
southwestern shell at 12.2–15.0 km s-1. Moreover, 13CO clumps at 3.8–9.4 km
s-1 spatially coincide well with radio filaments inside the shell. Such
spatial correspondence is naturally expected as a result of shock–cloud
interactions. According to magneto-hydrodynamical (MHD) simulations,
interactions between supernova shocks and clumpy clouds enhance turbulent
magnetic field up to $\sim$1 mG on the surface of shocked clouds, where the
synchrotron radio/X-ray radiation becomes brighter (Inoue et al., 2009, 2012;
Celli et al., 2019). This was further supported by several observations toward
the Galactic and Magellanic SNRs (Sano et al., 2013, 2017a, 2017b, 2019b,
2020b; Yamane et al., 2018; Kuriki et al., 2018), and hence it should not be
surprising that shock–cloud interactions with the magnetic field amplification
occurred in W49B as well.

However, we cannot rule out the possibility of shock-interaction with the 40
km s-1 and 60 km s-1 clouds from the spatial comparative studies alone. In
fact, molecular clouds at 38.0–40.0 km s-1 and 42.0–44.0 km s-1 show good
spatial correspondences with the southeastern half and southwestern shell of
the SNR, respectively (see Figure 3). The 63.0–65.0 km s-1 CO map also shows a
good anti-correlation with the radio shell, which is not inconsistent with the
picture of magnetic field amplification via the shock–cloud interaction.

####  4.1.2 Shock and Wind Induced Expanding Gas Motion

We argue that the cavity-like structures in the position-velocity diagrams at
the 10 km s-1 and 40 km s-1 clouds provide further supports for the shock
interactions (see Figures 4a–4d). Because such cavity-like structure toward a
SNR indicates an expanding gas motion, which is thought to be formed by a
combination of shock acceleration and strong gas winds from the progenitor
system: stellar winds from a high-mass progenitor or disk winds from a
progenitor system of post single-degenerate explosion. In the present study,
an expanding velocity $\Delta V$ is derived to be $\sim$6 km s-1 for the 10 km
s-1 cloud; and $\sim$3 km s-1 for the 40 km s-1 cloud. These values are
roughly consistent with other Galactic/Magellanic SNRs (e.g., Koo et al.,
1990; Koo & Heiles, 1991; Sano et al., 2017b, 2019c; Kuriki et al., 2018).

It is noteworthy that the two expanding cavities are independent because their
$\Delta V$ values are much smaller than the velocity difference of the 10 km
s-1 and 40 km s-1 clouds. Therefore, either expanding shell is located at the
same distance with W49B, and the forward shock has been impacted the wind-
cavity wall where the shock–cloud interaction occurred. The other expanding
shell is likely not associated with W49B. According to Sofue (2020), there are
many relics of fully evolved SNRs in the Galactic plane that cannot be
observed by radio continuum, optical, infrared, and X-rays. Because thermal
radiation from SNRs has shorter cooling time (below $\sim$10 kyr) than the
lifetime of giant molecular clouds ($\sim$10 Myr, e.g., Fukui et al., 1999;
Kawamura et al., 2009). The expanding gas motion of the 10 km s-1 or 40 km s-1
cloud is therefore likely one of such objects that happen to be located along
the line of sight.

####  4.1.3 Kinetic Temperature of Molecular Clouds

According to the LVG analysis in Section 3.5, the higher kinetic temperature
$T_{\mathrm{kin}}\sim$20–60 K of the shell clouds are seen in the 10 km s-1
cloud, suggesting that shock-heating likely occurred. Because the
$T_{\mathrm{kin}}$ values are roughly consistent with the previous studies of
shock-heated molecular clouds in the vicinity of middle-aged SNRs (e.g., Seta
et al., 1998; Gusdorf et al., 2012; Yoshiike et al., 2013; Anderl et al.,
2014; Dell’Ova et al., 2020). In addition, the presence of high temperature
dust components of $45\pm 4$ K and $151\pm 20$ K also support the shock-
heating scenario (Zhu et al., 2014). Moreover, the bright 24-micron emission
is detected in the southwestern shell where SW-edge of the 10 km s-1 cloud
shows the highest kinetic temperature of $\sim$60 K. It is noteworthy that
there are no other extra heating sources such as IRAS point sources or H II
regions toward the shell clouds (see also Figure 5).

By contrast, all shell clouds in the 40 km s-1 and 60 km s-1 components show
$T_{\mathrm{kin}}\sim$10 K, implying quiescent molecular clouds without any
extra-heating processes such as shock heating and stellar radiation.
Interestingly, Reference-SE in the 60 km s-1 cloud shows warmer temperature of
$\sim$20 K, despite the reference cloud is far from the SNR shell. A possible
scenario is that a part of the 60 km s-1 cloud is located at the tangent point
of the Galaxy, and hence the velocity crowding would accumulate diffuse gas
and increase the ambient gas temperature (e.g., Liu et al., 2019). In any
case, there is no shock-heated gas in both the 40 km s-1 and 60 km s-1 clouds.

Table 1: Results of LVG analysis at the 10 km s-1, 40 km s-1, and 60 km s-1
clouds

| $\mathrm{{}^{12}CO}$ |  | $\mathrm{{}^{13}CO}$ |  |
---|---|---|---|---|---
Name |  $J$ = 3–2 |  $J$ = 2–1 |  |  $J$ = 2–1 |  $n$(H2) |  0$T_{\mathrm{kin}}$0
| (K) | (K) |  | (K) |  0($\times 10^{3}$ cm-3)0 | (K)
(1) | (2) | (3) |  | (4) | (5) | (6)
10 km s-1 cloud000 |  |  |  |  |  |
SW-edge | 1.02 | 1.35 |  | 0.14 | $0.83^{+0.10}_{-0.05}$ | $60^{+47}_{-26}$
SW-shell | 3.12 | 3.97 |  | 1.12 | $2.45^{+0.50}_{-0.16}$ | ${18}^{{+7}}_{-3\phantom{0}}$
Inner-filament | 2.06 | 3.25 |  | 0.37 | $0.78^{+0.03}_{-0.02}$ | $23^{+9}_{-5\phantom{0}}$
NE-shell | 2.16 | 2.96 |  | 0.60 | $1.26^{+0.15}_{-0.09}$ | ${20}^{{+8}}_{-4\phantom{0}}$
Reference-SW | 3.02 | 4.03 |  | 1.43 | $1.86^{+0.14}_{-0.08}$ | ${15}^{{+3}}_{-2\phantom{0}}$
Reference-NW | 2.49 | 3.19 |  | 1.41 | $2.29^{+0.46}_{-0.20}$ | ${15}^{+6}_{{-3}\phantom{0}}$
40 km s-1 cloud00 |  |  |  |  |  |
SE-shell | 1.81 | 3.32 |  | 0.68 | $1.78^{+0.08}_{-0.00}$ | $10^{{+2}}_{{-1}\phantom{0}}$
SW-edge | 1.78 | 3.13 |  | 2.40 | $1.26^{+0.03}_{-0.03}$ | $\phantom{0}9^{+2}_{-1\phantom{0}}$
Inner-filament | 0.64 | 1.91 |  | 0.22 | $0.69^{+0.02}_{-0.00}$ | $\phantom{0}8^{+3}_{-2\phantom{0}}$
N-shell | 0.93 | 5.22 |  | 0.23 | $0.76^{+0.05}_{-0.02}$ | $13^{{+6}}_{{-3}\phantom{0}}$
Reference-SW | 1.65 | 2.83 |  | 0.77 | $1.32^{+0.09}_{-0.09}$ | $10^{+2}_{-2\phantom{0}}$
Reference-N | 1.05 | 2.40 |  | 0.42 | $0.93^{+0.03}_{-0.00}$ | $\phantom{0}9^{+2}_{-2\phantom{0}}$
60 km s-1 cloud00 |  |  |  |  |  |
SW-edge | 2.78 | 4.61 |  | 1.62 | $2.04^{+0.05}_{-0.04}$ | $\phantom{0}9^{+1}_{-1\phantom{0}}$
SW-shell | 2.09 | 3.66 |  | 1.19 | $1.48^{+0.03}_{-0.03}$ | $\phantom{0}9^{+1}_{-1\phantom{0}}$
Inner-filament | 5.22 | 7.93 |  | 2.82 | $1.66^{+0.08}_{-0.04}$ | $11^{{+1}}_{-1\phantom{0}}$
NW-shell | 3.98 | 6.13 |  | 1.31 | $1.51^{+0.04}_{-0.03}$ | ${14}^{+1}_{{-2}\phantom{0}}$
Reference-SW | 1.76 | 3.64 |  | 0.42 | $0.91^{+0.02}_{-0.00}$ | $12^{{+3}}_{{-1}\phantom{0}}$
Reference-SE | 5.20 | 7.15 |  | 1.33 | $1.17^{{+0.06}}_{{-0.02}}$ | $21^{+3}_{-2\phantom{0}}$

Note. — Col. (1): Region name for each cloud. Cols. (2)–(3): Radiation
temperature for each line emission derived by the least-squares fitting using
a single Gaussian function. Col. (4): Number density of molecular hydrogen.
Col. (5): Kinetic temperature.

####  4.1.4 Final Decision and Consistency with Previous Studies

In conclusion, we claim that the 10 km s-1 cloud is the one most likely
associated with W49B in terms of its spatial distribution, kinetics, and
physical conditions. This velocity is consistent not only with the line-
broadening measurements by Kilpatrick et al. (2016), but also with the latest
H I absorption measurement toward W49B by Ranasinghe & Leahy (2018). In this
case, the kinematic distance of W49B is slightly revised to $11.0\pm 0.4$ kpc
assuming the Galactic rotation curve model of Brand & Blitz (1993) and the
latest Galactic parameters of $R_{0}=7.92$ kpc and $\Theta_{0}=227$ km s-1
(VERA Collaboration et al., 2020). This value is also roughly consistent with
the previous distance to W49B of $11.3\pm 0.4$ kpc derived by H I absorption
(Ranasinghe & Leahy, 2018).

On the other hand, there is large gap in radial velocities between the 10 km
s-1 cloud and shocked H2 line emission at $\sim$64 km s-1 (Lee et al., 2020).
We argue that this inconsistency should not be a problem considering the
excitation condition for each line emission. In general, the CO line emission
at 2.6 mm (also known as 12CO $J$ = 1–0 transition line) can trace a bulk mass
of molecular cloud with low kinetic temperature of $\sim$10 K. On the other
hand, the supernova-shocked H2 line emission traces only a small portion of
molecular cloud which is highly excited into $\sim$2000–3000 K (e.g., Mouri,
1994; Lee et al., 2020). In W49B, the CO-traced molecular cloud mass is
$\sim$$2.7\times 10^{4}$ $M_{\sun}$ for the 10 km s-1 cloud, whereas the mass
of shocked H2 is only 14–550 $M_{\sun}$ (Keohane et al., 2007). This indicates
that the shocked H2 mass is only $\sim$2% of the CO-traced molecular cloud
mass at most. Considering the momentum conservation between the SNR shocks and
interacting clouds, the shocked H2 component is more easily accelerated than
the CO-traced molecular cloud. We therefore argue that the velocity
inconsistency between the CO cloud and shocked H2 component is caused by the
difference of their excitation conditions and masses.

In any case, the physical interaction of the 10 km s-1 cloud with W49B means
that the bulk mass of molecular clouds is concentrated in the northwestern
half of W49B, not in the southwest shell. This is consistent with the hydrogen
density maps derived by Zhou & Vink (2018), who found higher plasma densities
in the west of W49B by the X-ray spectral modeling. The inhomogeneous gas
distribution will significantly affect to understand the origins of
recombining plasma and gamma-rays from W49B. We will discuss later them the
latter Sections 4.3 and 4.4.

###  4.2 A Detailed Comparison with X-Rays

To reveal a physical relation between the 10 km s-1 cloud and X-ray radiation,
we here compare the CO distributions with Chandra X-ray images222Because the
shell boundary of X-rays is almost similar to that of radio continuum except
for the western half, we here only present a spatial comparison with the CO
map of 12.2–15.0 km s-1 which is bright in the western part of the shell..
Figure 8a shows the ALMA ACA 12CO($J$ = 2–1) integrated intensity overlaid
with the Chandra X-ray contours at the energy band of 2–7 keV. In the
integration velocity range of 12.2–15.0 km s-1, CO clouds are perfectly along
with the zigzag pattern of the western X-ray shell, indicating that the shock
ionization occurred. Note that the spatial correspondence is also seen in the
X-ray image at 4.2–5.5 keV, which are mostly free from the interstellar
absorption. We thus suggest that the shockwave was strongly decelerated and
deformed in the western shell along the dense clouds, whereas the eastern
shell was freely expanded with a smooth shape of the forward shock. This also
indicates that the shock velocity of the eastern shell is faster than that of
the western shell. Further proper motion studies might be able to reveal the
velocity difference in the east–west direction.

Figure 8: (a) Integrated intensity maps of ALMA ACA 12CO($J$ = 2–1) superposed
on the Chandra X-ray intensity contours in the energy band of 2–7 keV. The
integration velocity range of CO is from 12.2 to 15.0 km s-1. The contour
levels are 0.3, 0.4, 0.7, 1.2, 1.9, 2.8, 3.9, and $5.2\times 10^{-7}$ photon
cm-2 s-1. (b) Integrated intensity maps of ALMA ACA 13CO($J$ = 2–1) superposed
on the continuum-subtracted Fe He$\alpha$ emission. The integration velocity
range of CO is from 1.0 to 15.0 km s-1. The contour levels are 0.5, 0.8, 1.1,
1.4, 1.7, and $2.0\times 10^{-7}$ photon cm-2 s-1.

Figure 8b shows an overlay map of the 13CO($J$ = 2–1) intensity image and the
continuum-subtracted Fe He$\alpha$ emission in white contours. To compare the
Fe-rich ejecta with the total amount of dense clouds, we use 13CO with the
whole velocity range of 1.0–15.0 km s-1. Although the elongated structure of
Fe-rich ejecta is believed to be related to a bipolar/jet-driven Type Ib/Ic
explosion (Keohane et al., 2007; Lopez et al., 2013a), the Fe-rich ejecta is
mainly located on the void of dense molecular clouds. Moreover, almost Fe-rich
ejecta is surrounded by dense molecular clumps. We argue that this situation
is consistent with the supernova explosion inside a barrel-shaped cavity which
was proposed by Zhou & Vink (2018). The authors revealed that an enhancement
of cool plasma component along the Fe-rich ejecta (or the void of dense
clouds) was observed by a spatially resolved X-ray spectroscopy (see Figure 4
in Zhou & Vink 2018). Following the proposed scenario, the forward shock was
freely expanded in the low-density medium at the beginning, and then suddenly
encountered with the dense gaseous materials traced by 13CO line emission
and/or cool plasma component. Since the shock–cloud interaction generates
multiple reflected (or inward) shocks, the Fe-rich ejecta is efficiently
heated up at higher densities toward the center of the SNR (see also Sano et
al., 2019b). The X-Ray Imaging and Spectroscopy Mission (XRISM, XRISM Science
Team, 2020) will provide us with further understanding of shock-interactions
through a detailed spatial comparison between the X-ray derived ionic
properties and CO clouds.

Figure 9: (a) Integrated intensity maps of ALMA ACA 12CO($J$ = 2–1) for the 10
km s-1 cloud superposed on NuSTAR Fe He$\alpha$ flux contours (Yamaguchi et
al., 2018). The integration velocity range and contour levels are the same as
shown in Figures 2a and 6, respectively. The bashed 12 boxes indicate
$1^{\prime}\times 1^{\prime}$ regions used for the spatially resolved spectral
analysis in Yamaguchi et al. (2018) and deriving CO averaged integrated
intensities in Figure 9b. (b) Scatter plot between the electron temperature
$kT_{\mathrm{e}}$ (Yamaguchi et al., 2018) and peak integrated intensities of
12CO($J$ = 2–1) for box regions A1–3, B1–3, C1–3, and D1–3 as shown in Figure
9a. Error bars of CO and $kT_{\mathrm{e}}$ represent standard division of CO
integrated intensity and 1$\sigma$ confidence level for each box region. The
dashed line indicates the linear regression applying the least squares method.

###  4.3 Origin of the High-Temperature Recombining Plasma in W49B

It is a long-standing question how the recombining (overionized) plasma is
formed in SNRs since its discovery in 2002 (IC443, Kawasaki et al., 2002).
Subsequent detailed X-ray spectroscopic observations revealed that nearly 20
SNRs show the overionized state (e.g., W49B, Ozawa et al. 2009; G359.1$-$0.5,
Ohnishi et al. 2011; W28, Sawada & Koyama 2012; W44, Uchida et al. 2012;
G346.6$-$0.2, Yamauchi et al. 2013; 3C 391, Ergin et al. 2014; CTB 37A,
Yamauchi et al. 2014; G290.1$-$0.8, Kamitsukasa et al. 2015; LMC N49, Uchida
et al. 2015; Kes 17, Washino et al. 2016; G166.0$+$4.3, Matsumura et al.
2017a; 3C400.2, Ergin et al. 2017; LMC N132D, Bamba et al. 2018; HB21, Suzuki
et al. 2018; CTB1, Katsuragawa et al. 2018; Sagittarius A East, Ono et al.
2019; G189.6+3.3, Yamauchi et al. 2020, see also a review by Yamaguchi 2020).
However, the physical origin of recombining plasmas is still under debate.

Since the recombining plasma is characterized by higher ionization temperature
$kT_{\mathrm{i}}$ than the electron temperature $kT_{\mathrm{e}}$, rapid
electron cooling or increasing ionization state is needed to produce the
plasma state. Three scenarios have been proposed to explain the origin of
recombining plasmas in SNRs, called adiabatic cooling, thermal conduction, and
photoionization scenarios. In the adiabatic cooling (a.k.a. rarefaction)
scenario, rapid electron cooling occurs when the shockwaves breakout from a
dense ISM (e.g., CSM) into a much less dense medium (e.g., Itoh & Masai, 1989;
Masai, 1994; Yamaguchi et al., 2018). In the thermal conduction scenario, such
rapid electron cooling is caused by interactions between the shockwaves and
cold dense clouds through thermal conduction (e.g., Kawasaki et al., 2002;
Matsumura et al., 2017a, b; Okon et al., 2018, 2020). On the other hand, the
photoionization scenario proposes that an external X-ray radiation or low-
energy CRs) increase the ionization state via photoionization (e.g., Nakashima
et al., 2013; Ono et al., 2019; Hirayama et al., 2019). Because
photoionization can be seen in limited environments such as near the Galactic
center or a SNR with strong Fe I K$\alpha$ emission, the adiabatic cooling and
thermal conduction scenarios are thought to be the formation mechanisms of
recombining plasmas in most SNRs.

The origin of recombining plasma in W49B has been discussed in the past
decade. The thermal conduction scenario was initially proposed by Kawasaki et
al. (2005), whereas the adiabatic cooling scenario is more favored in the
subsequent studies (Miceli et al., 2010; Lopez et al., 2013a; Zhou et al.,
2011; Yamaguchi et al., 2018). Because the recombining plasma in W49B shows a
positive correlation between the ionization timescale $n_{\mathrm{e}}t$ and
$kT_{\mathrm{e}}$. Further, there is no correlation between the plasma
condition and ambient clouds traced by near infrared emission (Yamaguchi et
al., 2018). This trend is in contrast to what is observed in W44 (see also
Okon et al. 2020 and a review in Yamaguchi 2020). On the other hand, most
recent X-ray studies presented that the X-ray spectra from W49B are reproduced
by two ejecta components (low- and high-temperature plasma). The authors
proposed thermal conduction scenarios especially for the high-temperature
recombining plasma in W49B, considering the conduction timescale (Sun & Chen,
2020; Holland-Ashford et al., 2020). Note that Holland-Ashford et al. (2020)
argued that thermal conduction is a possible origin of recombining plasma in
the eastern regions of W49B because dense molecular clouds are thought to be
associated in the southwestern shell (Keohane et al., 2007; Zhu et al., 2014).
In this section, we argue that the origin of the high-temperature recombining
plasma in W49B can be understood as the thermal conduction scenario
considering the CO-traced interacting molecular clouds in W49B.

Table 2: Comparison of Physical Properties in Eleven Gamma-Ray SNRs

Name | Distance | Diameter | Age | $n_{\mathrm{p}}$ | $W_{\mathrm{p}}$ | References
---|---|---|---|---|---|---
| (kpc) | (pc) | (kyr) | (cm-3) | ($10^{49}$ erg) |
(1) | (2) | (3) | (4) | (5) | (6) | (7)
RX J1713.7$-3946$ | 1.0 | 18 | 1.6 | 130 |  00$0.16^{+0.07}_{-0.08}$ | Fukui et al. (2012)
RX J0852.0$-$4622 |  0z0.75a | 24 |  z1.7a | 100 |  00$0.07^{+0.02}_{-0.02}$ | Fukui et al. (2017)
RCW 86 | 2.5 | 30 | 1.8 |  075 |  00$0.11^{+0.01}_{-0.01}$ | Sano et al. (2019a)
HESS J1731$-$347 | 5.7 | 44 | 4.0 |  060 |  00$0.66^{+0.22}_{-0.22}$ | Fukuda et al. (2014)
G39.2$-$0.3 | 6.2 | 14 |  00zw$5.0^{+2.0}_{-2.0}$b | 400 | $3.2^{+1.1}_{-0.8}$ | de Oña Wilhelmi et al. (2020)
W49B | 11.00 | 16 |  00zw$6.0^{+1.0}_{-1.0}$c | 650 | $2.1^{+1.1}_{-0.6}$ | This work
Kes 79 | 5.5 | 16 |  0zw$8.3^{+0.5}_{-0.5}$ | 360 | 0.50zw | Kuriki et al. (2018)
W44 |  03.0d | 27 | 20.0e | 200 | 1.00zw | Yoshiike et al. (2013)
IC443 |  i1.5f | 20 |  0iiiii$25.0^{+5.0}_{-5.0}$g | 680 | 0.09ccz | Yoshiike et al. (2021)
LMC N132D | 50.00 | 25 |  00zw$2.5^{+0.2}_{-0.2}$h | $<2000$ | $>0.5$ | Sano et al. (2020a)
LMC N63A | 50.00 | 18 |  00zii$3.5^{+1.5}_{-1.5}$i | 190 | $0.9^{+0.5}_{-0.6}$ | Sano et al. (2019b)

Note. — Col. (1): Name of SNRs. Col. (2): Distance to SNRs in units of kpc.
Col. (3): Diameter of SNRs in units of pc. Col. (4): Age of SNRs in units of
kyr. Col. (5): Averaged number density of total interstellar protons
$n_{\mathrm{p}}$ in units of cm-3. Col. (6): Total energy of cosmic-ray
protons $W_{\mathrm{p}}$ in units of 1049 erg. Col. (7): References for CO/H I
derived $n_{\mathrm{p}}$ and $W_{\mathrm{p}}$ for each SNR. Other specific
references are also shown as follows: aKatsuda et al. (2008), bSu et al.
(2011), cZhou & Vink (2018), dCaswell et al. (1975), eWolszczan et al. (1991),
fWelsh & Sallmen (2003), gLee et al. (2008); Olbert et al. (2001), hLaw et al.
(2020), and iHughes et al. (1998).

Figure 9a shows the 12CO($J$ = 2–1) integrated intensity map of the 10 km s-1
cloud superposed on the NuSTAR Fe He$\alpha$ flux contours (Yamaguchi et al.,
2018). The twelve $1^{\prime}\times 1^{\prime}$ boxes indicate regions used
for spatially resolved spectral analysis using NuSTAR by Yamaguchi et al.
(2018). The authors fitted X-ray spectra above 3 keV using a single
temperature model, because the energy band is dominated by the high-
temperature plasma. We note that the CO integrated intensities are
significantly changed region to region, which shows an intensity gradient from
the southeast to the northwest as mentioned in Section 3.

Figure 9b shows a scatter plot between $kT_{\mathrm{e}}$ of the high-
temperature recombining plasma (Yamaguchi et al., 2018) and peak integrated
intensity of 12CO($J$ = 2–1) line emission for each box. We find a clear
negative correlation between the two. More precisely, $kT_{\mathrm{e}}$ values
in high-temperature plasma are increasing from the west (cloud rich) regions
to the east (cloud poor) regions. This is consistent with the thermal
conduction scenario: rapid electron cooling occurred in cold/dense cloud rich
regions. Note that this finding will not rule out the adiabatic cooling
scenario in W49B. In fact, the X-ray spectra from W49B are reproduced two
ejecta components: the low-temperature recombining plasma favors the adiabatic
cooling scenario whereas the high-temperature component is likely produced by
the thermal conduction (Sun & Chen, 2020; Holland-Ashford et al., 2020). In
other words, both the thermal conduction and adiabatic cooling processes
coexist in W49B.

Finally, we discuss the reason why our conclusion—the thermal conduction
origin of the high-temperature plasma—is different from some previous studies.
One of the most important issues is the previous evaluation of the ISM
interacting with W49B. Almost all previous studies used the shocked H2
distribution as the bulk mass of the ISM. However, as discussed in Section
4.1.4, the shocked H2 mass is only $\sim$2% of the CO-traced molecular cloud
mass. Because the shocked H2 map is bright in the southeast, most of
researchers believed that the southeast shell is interacting with dense
molecular clouds and the ISM mass of east is higher than that of west. Some
previous studies therefore concluded that the lower plasma temperature in west
was caused by the adiabatic cooling process (e.g., Holland-Ashford et al.,
2020; Yamaguchi et al., 2018). By contrast, it is noteworthy that Zhou & Vink
(2018) suggested that molecular cloud density is higher in the west than the
east, which is compatible with our ALMA results. We also note there are
different interpretations for the $n_{\mathrm{e}}t$ variation in the high-
temperature plasma. Yamaguchi et al. (2018) used $n_{\mathrm{e}}t$ as a proxy
for electron density, assuming a uniform time since heating and uniform
initial temperature. The former contrasts with the results of Holland-Ashford
et al. (2020) and Zhou & Vink (2018), who found higher recombination ages in
the east than the west. The positive correlation between $n_{\mathrm{e}}t$ and
$kT_{\mathrm{e}}$ in W49B may have to be reconsidered. In any case, we
emphasize that the proper evaluation of the ISM surrounding an SNR is
essential to understand the origin of recombining plasma correctly. Further
detailed comparative studies of CO based molecular cloud properties and X-ray
spectroscopic results are needed to better understand the origin of
recombining plasma in SNRs.

Figure 10: Scatter plot between the age of SNRs and the total energy of
cosmic-ray protons $W_{\mathrm{p}}$. The data points and references are
summarized in Table 2. The green solid line indicates the linear regression of
the double-logarithmic plot applying the least squares method.

###  4.4 Total Energy of Cosmic-Ray Protons

It is a hundred-year problem how CRs, mainly comprising relativistic protons,
are accelerated in interstellar space. SNRs are believed to be acceleration
sites for Galactic CRs below $\sim$3 PeV through the diffusive shock
acceleration (DSA, e.g., Bell, 1978; Blandford & Ostriker, 1978). A
conventional value of the total energy of CRs accelerated in an SNR is thought
to be $\sim$$10^{49}$–$10^{50}$ erg, corresponding to $\sim$1–10% of typical
kinematic energy released by a supernova ($10^{51}$ erg, e.g., Gabici, 2013;
Leahy et al., 2019). One of the foremost challenges is to validate these
predictions experimentally.

Gamma-ray and radio-line observations hold a key to understand the
acceleration of CRs in SNRs. Gamma-rays from SNRs are produced by two
different mechanisms: hadronic and leptonic processes (e.g., Aharonian et al.,
1994; Drury et al., 1994). For the hadronic process, CR proton–interstellar
proton interaction creates a neutral pion that quickly decays into two gamma-
ray photons (hadronic gamma-ray). For the leptonic scenario, a CR electron
energizes a low-energy photon into gamma-ray energy via inverse Compton
scattering, in addition to produce gamma-rays through non-thermal
Bremsstrahlung (leptonic gamma-ray). To confirm the acceleration of
relativistic protons, the main component of CRs, it is crucial to detect the
characteristic spectral feature of hadronic gamma-rays with a cut-off at a few
GeV known as pion-decay bump (e.g., Giuliani et al., 2011; Ackermann et al.,
2013). In addition, a good spatial correspondence between gamma-rays and
interstellar protons provides an alternative support for the CR proton
acceleration (e.g., Fukui et al., 2003, 2012, 2017; Aharonian et al., 2008;
Yoshiike et al., 2013; Sano et al., 2019a), because the hadronic gamma-ray
flux is proportional to the total energy of CR protons and the number density
of interstellar protons. Note that the interstellar protons are mainly either
neutral molecular and atomic hydrogen traced by CO and H I radio-lines,
respectively. Adopting the number density of targeted interstellar protons,
the total energy of CR protons was estimated to be
$W_{\mathrm{p}}\sim$1048–$10^{49}$ erg toward a dozen gamma-ray SNRs. However,
it is still under debate which parameters are important to understand the
variety of observed (or in-situ) $W_{\mathrm{p}}$ values. To better understand
the origins of CR protons and variety of $W_{\mathrm{p}}$, we need more
samples as well as detailed gamma-ray and radio-line studies for SNRs.

W49B is thought to be one of the CR proton accelerators because of the
detection of hadron-dominant gamma-rays with a pion-decay bump (H. E. S. S.
Collaboration et al., 2018). In fact, the best-fit position of GeV gamma-ray
detected by Fermi Large Area Telescope (LAT) is the edge of the SE shell,
where the bright radio continuum, shocked H2 emission, and the dense molecular
clouds are located. To obtain the total energy of CR protons in W49B, we first
estimate the number density of interstellar protons interacting with the SNR.
Using the equations (1) and (2), the averaged number density of interstellar
protons in molecular form is estimated to be $\sim$$650\pm 200$ cm-3 assuming
a shell radius of 8 pc and a thickness of 3 pc (e.g., Moffett & Reynolds,
1994). The error is derived as the typical uncertainty of the CO-to-H2
conversion factor of $\sim$30% (cf. Bolatto et al., 2013). Additionally, the
interstellar protons in atomic form are neglectable in W49B because the
derived column density of atomic hydrogen is significantly lower than that of
molecular hydrogen (see Brogan & Troland, 2001). The similar situation is also
seen in other middle-aged SNRs (e.g., W44, Yoshiike et al., 2013). We
therefore adopt the number density of interstellar protons $n_{\mathrm{p}}$ to
be $\sim$$650\pm 200$ cm-3.

According to H. E. S. S. Collaboration et al. (2018), the total energy of CR
protons $W_{\mathrm{p}}$ is written as

$\displaystyle W_{\mathrm{p}}\sim 2.0\mathrm{-}2.2\times
10^{49}(n_{\mathrm{p}}/650\;\mathrm{cm^{-3}})^{-1}\;(d/11\;\mathrm{kpc})^{-2}\;\;\mathrm{erg},\;\;\;\;\;$
(3)

where $d$ is the distance to the SNR. Adopting $n_{\mathrm{p}}=650$ cm-3 and
$d=11$ kpc, we than obtain $W_{\mathrm{p}}\sim$$2\times 10^{49}$ erg,
corresponding to $\sim$2% of the typical kinematic energy released by a
supernova explosion. Table 2 compares physical properties of eleven gamma-ray
SNRs including W49B. Here, all values of $n_{\mathrm{p}}$ and $W_{\mathrm{p}}$
were derived from CO/H I radio-line observations. We find that
$W_{\mathrm{p}}$ in W49B is roughly consistent with that in other gamma-ray
bright SNRs located in our Galaxy or the Large Magellanic Cloud (LMC). In
addition, it is noteworthy that young SNRs RX J1713.7$-$3946, RX
J0852.0$-$4622 (a.k.a. Vela Jr.), and RCW 86 as well as an evolved SNR IC 443
show the lowest values of $W_{\mathrm{p}}\sim$$10^{48}$ erg, while the others
hold higher values of $W_{\mathrm{p}}\sim$$10^{49}$ erg.

To better understand the trend, we plot $W_{\mathrm{p}}$ values as a function
of the age of SNRs. Figure 10 shows a scatter plot between the age of SNRs and
$W_{\mathrm{p}}$. We find a positive correlation between two parameters in the
SNRs with a young age less than $\sim$6000 yrs, suggesting that in-situ values
of $W_{\mathrm{p}}$ are strongly limited by short duration time of
acceleration also known as age-limited acceleration (cf. Ohira et al., 2010).
On the other hand, other SNRs with an older age more than $\sim$8000 yrs show
a steady decrease of $W_{\mathrm{p}}$ as SNRs get older. This trend could be
understood considering the energy dependent diffusion of CRs (e.g., Aharonian
& Atoyan, 1996; Gabici et al., 2007). In other words, in-situ values of
$W_{\mathrm{p}}$ have been decreased due to CR escape from the SNR. In fact,
hadron-dominant gamma-rays have been detected in nearby giant molecular clouds
of W44, suggesting that the molecular clouds are illuminated by CR protons
escaped from W44 (e.g., Uchiyama et al., 2012; Peron et al., 2020). The
authors suggested an actual value of $W_{\mathrm{p}}$ including escaped CRs is
$\sim$$10^{50}$ erg, corresponding to 10% of the typical kinematic energy
released by a supernova explosion. In any case, W49B shows one of the highest
in-situ values of $W_{\mathrm{p}}$ in the gamma-ray bright SNRs, which imply
that the escape (diffusion) of CRs is not significant at the moment. Further
gamma-ray observations using the Cherenkov Telescope Array (CTA) will unveil a
transition phase from the age-limited acceleration to escape dominant stage in
detail.

##  5 CONCLUSIONS

We summarize the primary conclusions as follows:

  1. 1.

New ALMA ACA CO($J$ = 2–1) observations at $\sim$$7^{\prime\prime}$ resolution
have revealed the spatial and kinematic distributions of three candidates of
interacting molecular clouds with the mixed-morphorogy SNR W49B, velocities of
which are $\sim$10 km s-1, $\sim$40 km s-1, and $\sim$60 km s-1. We found that
western molecular clouds at $\sim$10 km s-1 are obviously along with both the
radio continuum boundary and inside filaments as well as the deformed X-ray
shell, suggesting that shock-cloud interactions occurred. The 10 km s-1 cloud
also shows higher kinetic temperature of $\sim$20–60 K than the reference
clouds at 15 K, indicating that modest shock heating also occurred. The
presence of a wind-bubble with an expanding velocity of $\sim$6 km s-1
provides further evidence for the association of the 10 km s-1 cloud.

  2. 2.

The barrel-like structure of Fe-rich ejecta is mainly located on the void of
dense molecular clouds, where a cool plasma component is enhanced. We propose
a possible scenario that the barrel-like structure of Fe-rich ejecta was
formed not only by the asymmetric supernova explosion, but also by
interactions with dense molecular clouds. A supernova explosion occurred
within the cylinder-like gaseous medium and then Fe-rich ejecta was
efficiently heated-up at higher densities by multiple-reflected shocks formed
by shock-cloud interactions.

  3. 3.

The electron temperature $kT_{\mathrm{e}}$ of recombining plasma from Fe
He$\alpha$ shows a negative correlation with the peak integrated intensity of
CO line emission in the 10 km s-1 cloud. More precisely, $kT_{\mathrm{e}}$
values in high-temperature recombining plasma are increasing from the west
(cloud rich) regions to the east (cloud poor) regions, suggesting the thermal
conduction origin. Note that this finding does not rule out the adiabatic
cooling scenario in the low-temperature recombining plasma in W49B which was
previously discussed (Sun & Chen, 2020; Holland-Ashford et al., 2020).

  4. 4.

The total energy of CR protons $W_{\mathrm{p}}$ is estimated to be
$\sim$$2\times 10^{49}$ erg, which is one of the highest values in gamma-ray
bright SNRs. We found that in-situ values of $W_{\mathrm{p}}$ in gamma-ray
SNRs increase with age for the young group (with the age less than $\sim$6000
yr). On the other hand, other older SNRs show a steady decrease of
$W_{\mathrm{p}}$ as SNRs get older due to the escapes/diffusion effect of CRs.
We frame a hypothesis that W49B is undergoing an age-limited acceleration
without a significant escape or diffusion of CRs from the SNR.

## Acknowledgements

We are grateful Hiroya Yamaguchi for providing us with the NuSTAR data points
used in this paper. This paper makes use of the following ALMA data:
ADS/JAO.ALMA#2018.1.01780.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. The scientific results reported in this article are based on data
obtained from the Chandra Data Archive (Obs IDs: 117, 13440, and 13441). This
research has made use of software provided by the Chandra X-ray Center (CXC)
in the application packages CIAO (v 4.12). This work was supported by JSPS
KAKENHI Grant Numbers JP19H05075 (H. Sano), and JP21H01136 (H. Sano). K.
Tokuda was supported by NAOJ ALMA Scientific Research Grant Number of
2016-03B. We are also grateful to the anonymous referee for useful comments
which helped the authors to improve the paper significantly.

## References

  * Abdo et al. (2010) Abdo, A. A., Ackermann, M., Ajello, M., et al. 2010, ApJ, 722, 1303. doi:10.1088/0004-637X/722/2/1303
  * Ackermann et al. (2013) Ackermann, M., Ajello, M., Allafort, A., et al. 2013, Science, 339, 807. doi:10.1126/science.1231160
  * Aharonian et al. (1994) Aharonian, F. A., Drury, L. O., & Voelk, H. J. 1994, A&A, 285, 645
  * Aharonian & Atoyan (1996) Aharonian, F. A. & Atoyan, A. M. 1996, A&A, 309, 917
  * Aharonian et al. (2008) Aharonian, F., Akhperjanian, A. G., Bazer-Bachi, A. R., et al. 2008, A&A, 481, 401. doi:10.1051/0004-6361:20077765
  * Anderl et al. (2014) Anderl, S., Gusdorf, A., & Güsten, R. 2014, A&A, 569, A81. doi:10.1051/0004-6361/201423561
  * Bamba et al. (2016) Bamba, A., Sawada, M., Nakano, Y., et al. 2016, PASJ, 68, S5. doi:10.1093/pasj/psv096
  * Bamba et al. (2018) Bamba, A., Ohira, Y., Yamazaki, R., et al. 2018, ApJ, 854, 71. doi:10.3847/1538-4357/aaa5a0
  * Bear & Soker (2017) Bear, E. & Soker, N. 2017, MNRAS, 468, 140. doi:10.1093/mnras/stx431
  * Bieging & Peters (2011) Bieging, J. H. & Peters, W. L. 2011, ApJS, 196, 18. doi:10.1088/0067-0049/196/2/18
  * Beichman et al. (1988) Beichman, C. A., Neugebauer, G., Habing, H. J., et al. 1988, Infrared astronomical satellite (IRAS) catalogs and atlases. Volume 1: Explanatory supplement, 1
  * Bell (1978) Bell, A. R. 1978, MNRAS, 182, 147. doi:10.1093/mnras/182.2.147
  * Bertsch et al. (1993) Bertsch, D. L., Dame, T. M., Fichtel, C. E., et al. 1993, ApJ, 416, 587. doi:10.1086/173261
  * Blandford & Ostriker (1978) Blandford, R. D. & Ostriker, J. P. 1978, ApJ, 221, L29. doi:10.1086/182658
  * Blake et al. (1987) Blake, G. A., Sutton, E. C., Masson, C. R., et al. 1987, ApJ, 315, 621. doi:10.1086/165165
  * Bolatto et al. (2013) Bolatto, A. D., Wolfire, M., & Leroy, A. K. 2013, ARA&A, 51, 207. doi:10.1146/annurev-astro-082812-140944
  * Brand & Blitz (1993) Brand, J. & Blitz, L. 1993, A&A, 275, 67
  * Brogan & Troland (2001) Brogan, C. L. & Troland, T. H. 2001, ApJ, 550, 799. doi:10.1086/319787
  * Caswell et al. (1975) Caswell, J. L., Murray, J. D., Roger, R. S., et al. 1975, A&A, 45, 239
  * Celli et al. (2019) Celli, S., Morlino, G., Gabici, S., et al. 2019, MNRAS, 487, 3199. doi:10.1093/mnras/stz1425
  * Dell’Ova et al. (2020) Dell’Ova, P., Gusdorf, A., Gerin, M., et al. 2020, A&A, 644, A64. doi:10.1051/0004-6361/202038339
  * Dempsey et al. (2013) Dempsey, J. T., Thomas, H. S., & Currie, M. J. 2013, ApJS, 209, 8. doi:10.1088/0067-0049/209/1/8
  * de Oña Wilhelmi et al. (2020) de Oña Wilhelmi, E., Sushch, I., Brose, R., et al. 2020, MNRAS, 497, 3581. doi:10.1093/mnras/staa2045
  * Drury et al. (1994) Drury, L. O., Aharonian, F. A., & Voelk, H. J. 1994, A&A, 287, 959
  * Ergin et al. (2014) Ergin, T., Sezer, A., Saha, L., et al. 2014, ApJ, 790, 65. doi:10.1088/0004-637X/790/1/65
  * Ergin et al. (2017) Ergin, T., Sezer, A., Sano, H., et al. 2017, ApJ, 842, 22. doi:10.3847/1538-4357/aa72ee
  * Fukuda et al. (2014) Fukuda, T., Yoshiike, S., Sano, H., et al. 2014, ApJ, 788, 94. doi:10.1088/0004-637X/788/1/94
  * Fukui et al. (1999) Fukui, Y., Mizuno, N., Yamaguchi, R., et al. 1999, PASJ, 51, 745. doi:10.1093/pasj/51.6.745
  * Fukui et al. (2003) Fukui, Y., Moriguchi, Y., Tamura, K., et al. 2003, PASJ, 55, L61. doi:10.1093/pasj/55.5.L61
  * Fukui et al. (2012) Fukui, Y., Sano, H., Sato, J., et al. 2012, ApJ, 746, 82. doi:10.1088/0004-637X/746/1/82
  * Fukui et al. (2017) Fukui, Y., Sano, H., Sato, J., et al. 2017, ApJ, 850, 71. doi:10.3847/1538-4357/aa9219
  * Fruscione et al. (2006) Fruscione, A., McDowell, J. C., Allen, G. E., et al. 2006, Proc. SPIE, 6270, 62701V. doi:10.1117/12.671760
  * Gabici et al. (2007) Gabici, S., Aharonian, F. A., & Blasi, P. 2007, Ap&SS, 309, 365. doi:10.1007/s10509-007-9427-6
  * Gabici (2013) Gabici, S. 2013, Cosmic Rays in Star-Forming Environments, 34, 221. doi:10.1007/978-3-642-35410-6_16
  * Giuliani et al. (2011) Giuliani, A., Cardillo, M., Tavani, M., et al. 2011, ApJ, 742, L30. doi:10.1088/2041-8205/742/2/L30
  * Goldreich & Kwan (1974) Goldreich, P. & Kwan, J. 1974, ApJ, 189, 441. doi:10.1086/152821
  * Graessle et al. (2007) Graessle, D. E., Evans, I. N., Glotfelty, K., et al. 2007, Chandra Newsletter, Vol. 14, p.33, 14
  * Gusdorf et al. (2012) Gusdorf, A., Anderl, S., Güsten, R., et al. 2012, A&A, 542, L19. doi:10.1051/0004-6361/201218907
  * H. E. S. S. Collaboration et al. (2018) H. E. S. S. Collaboration, Abdalla, H., Abramowski, A., et al. 2018, A&A, 612, A5. doi:10.1051/0004-6361/201527843
  * Helfand et al. (2006) Helfand, D. J., Becker, R. H., White, R. L., et al. 2006, AJ, 131, 2525. doi:10.1086/503253
  * Hirayama et al. (2019) Hirayama, A., Yamauchi, S., Nobukawa, K. K., et al. 2019, PASJ, 71, 37. doi:10.1093/pasj/psy153
  * Holland-Ashford et al. (2020) Holland-Ashford, T., Lopez, L. A., & Auchettl, K. 2020, ApJ, 903, 108. doi:10.3847/1538-4357/abb808
  * Hughes et al. (1998) Hughes, J. P., Hayashi, I., & Koyama, K. 1998, ApJ, 505, 732. doi:10.1086/306202
  * Inoue et al. (2009) Inoue, T., Yamazaki, R., & Inutsuka, S.-. ichiro . 2009, ApJ, 695, 825. doi:10.1088/0004-637X/695/2/825
  * Inoue et al. (2012) Inoue, T., Yamazaki, R., Inutsuka, S.-. ichiro ., et al. 2012, ApJ, 744, 71. doi:10.1088/0004-637X/744/1/71
  * Itoh & Masai (1989) Itoh, H. & Masai, K. 1989, MNRAS, 236, 885. doi:10.1093/mnras/236.4.885
  * Kamitsukasa et al. (2015) Kamitsukasa, F., Koyama, K., Uchida, H., et al. 2015, PASJ, 67, 16. doi:10.1093/pasj/psu149
  * Katsuragawa et al. (2018) Katsuragawa, M., Nakashima, S., Matsumura, H., et al. 2018, PASJ, 70, 110. doi:10.1093/pasj/psy114
  * Kawasaki et al. (2002) Kawasaki, M. T., Ozaki, M., Nagase, F., et al. 2002, ApJ, 572, 897. doi:10.1086/340383
  * Kawasaki et al. (2005) Kawasaki, M., Ozaki, M., Nagase, F., et al. 2005, ApJ, 631, 935. doi:10.1086/432591
  * Kepley et al. (2020) Kepley, A. A., Tsutsumi, T., Brogan, C. L., et al. 2020, PASP, 132, 024505. doi:10.1088/1538-3873/ab5e14
  * Katsuda et al. (2008) Katsuda, S., Tsunemi, H., & Mori, K. 2008, ApJ, 678, L35. doi:10.1086/588499
  * Kawamura et al. (2009) Kawamura, A., Mizuno, Y., Minamidani, T., et al. 2009, ApJS, 184, 1. doi:10.1088/0067-0049/184/1/1
  * Keohane et al. (2007) Keohane, J. W., Reach, W. T., Rho, J., et al. 2007, ApJ, 654, 938. doi:10.1086/509311
  * Kerr & Lynden-Bell (1986) Kerr, F. J. & Lynden-Bell, D. 1986, MNRAS, 221, 1023. doi:10.1093/mnras/221.4.1023
  * Kilpatrick et al. (2016) Kilpatrick, C. D., Bieging, J. H., & Rieke, G. H. 2016, ApJ, 816, 1. doi:10.3847/0004-637X/816/1/1
  * Koo et al. (1990) Koo, B.-C., Reach, W. T., Heiles, C., et al. 1990, ApJ, 364, 178. doi:10.1086/169400
  * Koo & Heiles (1991) Koo, B.-C. & Heiles, C. 1991, ApJ, 382, 204. doi:10.1086/170709
  * Koo et al. (2016) Koo, B.-C., Lee, J.-J., Jeong, I.-G., et al. 2016, ApJ, 821, 20. doi:10.3847/0004-637X/821/1/20
  * Kuriki et al. (2018) Kuriki, M., Sano, H., Kuno, N., et al. 2018, ApJ, 864, 161. doi:10.3847/1538-4357/aad7be
  * Langer & Penzias (1990) Langer, W. D. & Penzias, A. A. 1990, ApJ, 357, 477. doi:10.1086/168935
  * Law et al. (2020) Law, C. J., Milisavljevic, D., Patnaude, D. J., et al. 2020, ApJ, 894, 73. doi:10.3847/1538-4357/ab873a
  * Leahy et al. (2019) Leahy, D., Wang, Y., Lawton, B., et al. 2019, AJ, 158, 149. doi:10.3847/1538-3881/ab3d2c
  * Lee et al. (2008) Lee, J.-J., Koo, B.-C., Yun, M. S., et al. 2008, AJ, 135, 796. doi:10.1088/0004-6256/135/3/796
  * Lee et al. (2020) Lee, Y.-H., Koo, B.-C., & Lee, J.-J. 2020, AJ, 160, 263. doi:10.3847/1538-3881/abc00e
  * Liu et al. (2019) Liu, B., Anderson, L. D., McIntyre, T., et al. 2019, ApJS, 240, 14. doi:10.3847/1538-4365/aaef8e
  * Lockhart & Goss (1978) Lockhart, I. A. & Goss, W. M. 1978, A&A, 67, 355
  * Lopez et al. (2009a) Lopez, L. A., Ramirez-Ruiz, E., Pooley, D. A., et al. 2009a, ApJ, 691, 875. doi:10.1088/0004-637X/691/1/875
  * Lopez et al. (2009b) Lopez, L. A., Ramirez-Ruiz, E., Badenes, C., et al. 2009b, ApJ, 706, L106. doi:10.1088/0004-637X/706/1/L106
  * Lopez et al. (2011) Lopez, L. A., Ramirez-Ruiz, E., Huppenkothen, D., et al. 2011, ApJ, 732, 114. doi:10.1088/0004-637X/732/2/114
  * Lopez et al. (2013a) Lopez, L. A., Pearson, S., Ramirez-Ruiz, E., et al. 2013a, ApJ, 777, 145. doi:10.1088/0004-637X/777/2/145
  * Lopez et al. (2013b) Lopez, L. A., Ramirez-Ruiz, E., Castro, D., et al. 2013b, ApJ, 764, 50. doi:10.1088/0004-637X/764/1/50
  * Masai (1994) Masai, K. 1994, ApJ, 437, 770. doi:10.1086/175037
  * Matsumura et al. (2017a) Matsumura, H., Uchida, H., Tanaka, T., et al. 2017a, PASJ, 69, 30. doi:10.1093/pasj/psx001
  * Matsumura et al. (2017b) Matsumura, H., Tanaka, T., Uchida, H., et al. 2017b, ApJ, 851, 73. doi:10.3847/1538-4357/aa9bdf
  * McMullin et al. (2007) McMullin, J. P., Waters, B., Schiebel, D., et al. 2007, Astronomical Data Analysis Software and Systems XVI, 376, 127
  * Miceli et al. (2010) Miceli, M., Bocchino, F., Decourchelle, A., et al. 2010, A&A, 514, L2. doi:10.1051/0004-6361/200913713
  * Moffett & Reynolds (1994) Moffett, D. A. & Reynolds, S. P. 1994, ApJ, 437, 705. doi:10.1086/175033
  * Mouri (1994) Mouri, H. 1994, ApJ, 427, 777. doi:10.1086/174184
  * Nakashima et al. (2013) Nakashima, S., Nobukawa, M., Uchida, H., et al. 2013, ApJ, 773, 20. doi:10.1088/0004-637X/773/1/20
  * Ohira et al. (2010) Ohira, Y., Murase, K., & Yamazaki, R. 2010, A&A, 513, A17. doi:10.1051/0004-6361/200913495
  * Ohnishi et al. (2011) Ohnishi, T., Koyama, K., Tsuru, T. G., et al. 2011, PASJ, 63, 527. doi:10.1093/pasj/63.3.527
  * Okon et al. (2018) Okon, H., Uchida, H., Tanaka, T., et al. 2018, PASJ, 70, 35. doi:10.1093/pasj/psy022
  * Okon et al. (2020) Okon, H., Tanaka, T., Uchida, H., et al. 2020, ApJ, 890, 62. doi:10.3847/1538-4357/ab6987
  * Olbert et al. (2001) Olbert, C. M., Clearfield, C. R., Williams, N. E., et al. 2001, ApJ, 554, L205. doi:10.1086/321708
  * Ono et al. (2019) Ono, A., Uchiyama, H., Yamauchi, S., et al. 2019, PASJ, 71, 52. doi:10.1093/pasj/psz025
  * Ozawa et al. (2009) Ozawa, M., Koyama, K., Yamaguchi, H., et al. 2009, ApJ, 706, L71. doi:10.1088/0004-637X/706/1/L71
  * Peron et al. (2020) Peron, G., Aharonian, F., Casanova, S., et al. 2020, ApJ, 896, L23. doi:10.3847/2041-8213/ab93d1
  * Ranasinghe & Leahy (2018) Ranasinghe, S. & Leahy, D. A. 2018, AJ, 155, 204. doi:10.3847/1538-3881/aab9be
  * Rho & Petre (1998) Rho, J. & Petre, R. 1998, ApJ, 503, L167. doi:10.1086/311538
  * Sano et al. (2013) Sano, H., Tanaka, T., Torii, K., et al. 2013, ApJ, 778, 59. doi:10.1088/0004-637X/778/1/59
  * Sano et al. (2017a) Sano, H., Yamane, Y., Voisin, F., et al. 2017a, ApJ, 843, 61. doi:10.3847/1538-4357/aa73e0
  * Sano et al. (2017b) Sano, H., Reynoso, E. M., Mitsuishi, I., et al. 2017b, Journal of High Energy Astrophysics, 15, 1. doi:10.1016/j.jheap.2017.04.002
  * Sano et al. (2019a) Sano, H., Rowell, G., Reynoso, E. M., et al. 2019a, ApJ, 876, 37. doi:10.3847/1538-4357/ab108f
  * Sano et al. (2019b) Sano, H., Matsumura, H., Nagaya, T., et al. 2019b, ApJ, 873, 40. doi:10.3847/1538-4357/ab02fd
  * Sano et al. (2019c) Sano, H., Matsumura, H., Yamane, Y., et al. 2019c, ApJ, 881, 85. doi:10.3847/1538-4357/ab2ade
  * Sano et al. (2020a) Sano, H., Plucinsky, P. P., Bamba, A., et al. 2020a, ApJ, 902, 53. doi:10.3847/1538-4357/abb469
  * Sano et al. (2020b) Sano, H., Inoue, T., Tokuda, K., et al. 2020b, ApJ, 904, L24. doi:10.3847/2041-8213/abc884
  * Sano & Fukui (2021) Sano, H. & Fukui, Y. 2021, arXiv:2106.00708
  * Sawada & Koyama (2012) Sawada, M. & Koyama, K. 2012, PASJ, 64, 81. doi:10.1093/pasj/64.4.81
  * Scoville & Solomon (1974) Scoville, N. Z. & Solomon, P. M. 1974, ApJ, 187, L67. doi:10.1086/181398
  * Seta et al. (1998) Seta, M., Hasegawa, T., Dame, T. M., et al. 1998, ApJ, 505, 286. doi:10.1086/306141
  * Siegel et al. (2020) Siegel, J., Dwarkadas, V. V., Frank, K. A., et al. 2020, ApJ, 904, 175. doi:10.3847/1538-4357/abbfa9
  * Shimizu et al. (2012) Shimizu, T., Masai, K., & Koyama, K. 2012, PASJ, 64, 24. doi:10.1093/pasj/64.2.24
  * Slavin et al. (2017) Slavin, J. D., Smith, R. K., Foster, A., et al. 2017, ApJ, 846, 77. doi:10.3847/1538-4357/aa8552
  * Sofue et al. (2020) Sofue, Y., Kohno, M., & Umemoto, T. 2020, arXiv:2011.11916
  * Sofue (2020) Sofue, Y. 2020, PASJ, 72, L11. doi:10.1093/pasj/psaa102
  * Su et al. (2011) Su, Y., Chen, Y., Yang, J., et al. 2011, ApJ, 727, 43. doi:10.1088/0004-637X/727/1/43
  * Sun & Chen (2020) Sun, L. & Chen, Y. 2020, ApJ, 893, 90. doi:10.3847/1538-4357/ab8001
  * Suzuki et al. (2018) Suzuki, H., Bamba, A., Nakazawa, K., et al. 2018, PASJ, 70, 75. doi:10.1093/pasj/psy069
  * Tanaka et al. (2018) Tanaka, T., Yamaguchi, H., Wik, D. R., et al. 2018, ApJ, 866, L26. doi:10.3847/2041-8213/aae709
  * Uchida et al. (2012) Uchida, H., Koyama, K., Yamaguchi, H., et al. 2012, PASJ, 64, 141. doi:10.1093/pasj/64.6.141
  * Uchida et al. (2015) Uchida, H., Koyama, K., & Yamaguchi, H. 2015, ApJ, 808, 77. doi:10.1088/0004-637X/808/1/77
  * Uchiyama et al. (2012) Uchiyama, Y., Funk, S., Katagiri, H., et al. 2012, ApJ, 749, L35. doi:10.1088/2041-8205/749/2/L35
  * Umemoto et al. (2017) Umemoto, T., Minamidani, T., Kuno, N., et al. 2017, PASJ, 69, 78. doi:10.1093/pasj/psx061
  * Urquhart et al. (2009) Urquhart, J. S., Hoare, M. G., Purcell, C. R., et al. 2009, A&A, 501, 539. doi:10.1051/0004-6361/200912108
  * VERA Collaboration et al. (2020) VERA Collaboration, Hirota, T., Nagayama, T., et al. 2020, PASJ, 72, 50. doi:10.1093/pasj/psaa018
  * Vink (2012) Vink, J. 2012, A&A Rev., 20, 49. doi:10.1007/s00159-011-0049-1
  * Washino et al. (2016) Washino, R., Uchida, H., Nobukawa, M., et al. 2016, PASJ, 68, S4. doi:10.1093/pasj/psv095
  * Welsh & Sallmen (2003) Welsh, B. Y. & Sallmen, S. 2003, A&A, 408, 545. doi:10.1051/0004-6361:20030908
  * Wolszczan et al. (1991) Wolszczan, A., Cordes, J. M., & Dewey, R. J. 1991, ApJ, 372, L99. doi:10.1086/186033
  * XRISM Science Team (2020) XRISM Science Team 2020, arXiv:2003.04962
  * Yamaguchi et al. (2014) Yamaguchi, H., Badenes, C., Petre, R., et al. 2014, ApJ, 785, L27. doi:10.1088/2041-8205/785/2/L27
  * Yamaguchi et al. (2018) Yamaguchi, H., Tanaka, T., Wik, D. R., et al. 2018, ApJ, 868, L35. doi:10.3847/2041-8213/aaf055
  * Yamaguchi (2020) Yamaguchi, H. 2020, Astronomische Nachrichten, 341, 150. doi:10.1002/asna.202023771
  * Yamauchi et al. (2014) Yamauchi, S., Minami, S., Ota, N., et al. 2014, PASJ, 66, 2. doi:10.1093/pasj/pst004
  * Yamane et al. (2018) Yamane, Y., Sano, H., van Loon, J. T., et al. 2018, ApJ, 863, 55. doi:10.3847/1538-4357/aacfff
  * Yamauchi et al. (2020) Yamauchi, S., Oya, M., Nobukawa, K. K., et al. 2020, PASJ, 72, 81. doi:10.1093/pasj/psaa070
  * Yamauchi et al. (2013) Yamauchi, S., Nobukawa, M., Koyama, K., et al. 2013, PASJ, 65, 6. doi:10.1093/pasj/65.1.6
  * Yang et al. (2009) Yang, X. J., Tsunemi, H., Lu, F. J., et al. 2009, ApJ, 692, 894. doi:10.1088/0004-637X/692/1/894
  * Yoshiike et al. (2013) Yoshiike, S., Fukuda, T., Sano, H., et al. 2013, ApJ, 768, 179. doi:10.1088/0004-637X/768/2/179
  * Yoshiike et al. (2021) Yoshiike, S., Sano, H., Fukuda, T., et al. 2021, to be submitted
  * Yusef-Zadeh et al. (2003) Yusef-Zadeh, F., Wardle, M., Rho, J., et al. 2003, ApJ, 585, 319. doi:10.1086/345932
  * Zhang et al. (2019) Zhang, G.-Y., Slavin, J. D., Foster, A., et al. 2019, ApJ, 875, 81. doi:10.3847/1538-4357/ab0f9a
  * Zhou et al. (2011) Zhou, X., Miceli, M., Bocchino, F., et al. 2011, MNRAS, 415, 244. doi:10.1111/j.1365-2966.2011.18695.x
  * Zhou & Vink (2018) Zhou, P. & Vink, J. 2018, A&A, 615, A150. doi:10.1051/0004-6361/201731583
  * Zhu et al. (2014) Zhu, H., Tian, W. W., & Zuo, P. 2014, ApJ, 793, 95. doi:10.1088/0004-637X/793/2/95