# Mapping and characterizing magnetic fields in the Rho Ophiuchus-A molecular cloud with SOFIA/HAWC+

Ngân Lê, Le Ngoc Tram, Agata Karska, Thiem Hoang, Pham Ngoc Diep, Michal Hanasz, Nguyen Bich Ngoc, Nguyen Thi Phuong, Karl M. Menten, Friedrich Wyrowski, et al.

<figure><img src="image_1.png" /></figure>

To cite this version:

Ngân Lê, Le Ngoc Tram, Agata Karska, Thiem Hoang, Pham Ngoc Diep, et al. Mapping and characterizing magnetic fields in the Rho Ophiuchus-A molecular cloud with SOFIA/HAWC+. Astronomy & Astrophysics - A&A, 2024, 690, 10.1051/0004-6361/202348008. insu-04838848

HAL Id: insu-04838848

https://insu.hal.science/insu-04838848v1

Submitted on 15 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Mapping and characterizing magnetic fields in the Rho Ophiuchus-A molecular cloud with SOFIA/HAWC+

Ngân Lê<fnref n="1" />,<fnref n="2" />,* , Le Ngoc Tram<fnref n="3" />, , Agata Karska<fnref n="1" />,<fnref n="3" />, , Thiem Hoang<fnref n="4" />,<fnref n="5" />, , Pham Ngoc Diep<fnref n="6" />,<fnref n="7" />,<fnref n="8" />, , Michał Hanasz<fnref n="1" />, , Nguyen Bich Ngoc<fnref n="6" />,<fnref n="7" />, , Nguyen Thi Phuong<fnref n="6" />,<fnref n="2" />, , Karl M. Menten<fnref n="3" />, , Friedrich Wyrowski<fnref n="3" />, , Dieu D. Nguyen<fnref n="8" />, , Thuong Duc Hoang<fnref n="9" />, , and Nguyen Minh Khang<fnref n="10" />

Received 18 September 2023 / Accepted 30 August 2024

## ABSTRACT

**Context.** Together with gravity, turbulence, and stellar feedback, magnetic fields (B-fields) are thought to play a critical role in the evolution of molecular clouds and star formation processes. The polarization of thermal dust emission is a popular tracer of B-fields in star-forming regions.

**Aims.** We aim to map the morphology and measure the strength of B-fields of the nearby molecular cloud, rho Ophiuchus-A ($\rho$ Oph-A), to understand the role of B-fields in regulating star formation and in shaping the cloud.

**Methods.** We analyzed the far-infrared (FIR) polarization of thermal dust emission observed by SOFIA/HAWC+ at 89 and 154 µm toward the densest part of $\rho$ Oph-A, which is irradiated by the nearby B3/4 star, Oph-S1. These FIR polarimetric maps cover an area of ~4.5' × 4.5' (corresponding to 0.18 × 0.18 pc<fnref n="2" />) with an angular resolution of 7.8" and 13.6", respectively.

**Results.** The $\rho$ Oph-A cloud exhibits well-ordered B-fields with magnetic orientations that are mainly perpendicular to the ridge of the cloud toward the densest region. We obtained a map of B-field strengths in the range of 0.2-2.5 mG, using the Davis-Chandrasekhar-Fermi (DCF) method. The B-fields are strongest at the densest part of the cloud, which is associated with the starless core SM1, and then decrease toward the outskirts of the cloud. By calculating the map of the mass-to-flux ratio, Alfvén Mach number, and plasma $\beta$ parameter in $\rho$ Oph-A, we find that the cloud is predominantly magnetically sub-critical, sub-Alfvénic, which implies that the cloud is supported by strong B-fields that dominate over gravity, turbulence, and thermal gas energy. The measured B-field strengths at the two densest subsregions using other methods that account for the compressible mode are relatively lower than that measured with the DCF method. However, these results do not significantly change our conclusions on the roles of B-fields relative to gravity and turbulence on star formation. Our virial analysis suggests that the cloud is gravitationally unbound, which is consistent with the previous detection of numerous starless cores in the cloud. By comparing the magnetic pressure with the radiation pressure from the Oph-S1 star, we find that B-fields are sufficiently strong to support the cloud against radiative feedback and to regulate the shape of the cloud.

**Key words.** stars: formation – ISM: clouds – ISM: magnetic fields – ISM: individual objects: Rho Oph-A

## 1. Introduction

Stars like our Sun form in the cold and dense parts of molecular clouds, which consist of filaments, clumps, and dense cores (Shu et al. 1987; Lada & Lada 2003; Kennicutt & Evans 2012). The formation of stars is regulated not only by gravity, but also by magnetic fields (B-fields), turbulence, and radiative and mechanical feedback from nearby stars (McKee & Ostriker 2007; Hennebelle & Falgarone 2012; Hennebelle & Inutsuka 2019).

Turbulence counteracts gravity and prevents the collapse of the dense core, but it can also generate shocks that compress the gas and trigger the formation of new stars.

B-fields are thought to play a critical role in the star formation process by supporting the material in the cloud against self-gravitational collapse (see, e.g., Mouschovias 2001; Crutcher 2012; Pattle & Fissel 2019; Pattle et al. 2023). The collapse can occur when the dense core is magnetically super-critical, with the ratio of mass to magnetic flux (known as mass-to-flux ratio) greater than 1; that is, when B-fields are not strong

<footnote n="1">¹ Institute of Astronomy, Faculty of Physics, Astronomy and Informatics, Nicolaus Copernicus University, Grudziładzka 5, 87-100 Toruń, Poland</footnote>

<footnote n="2">² Institute For Interdisciplinary Research in Science and Education (IFIRSE), ICISE, 07 Science Avenue, Ghenh Rang Ward, 55121 Quy Nhon City, Binh Dinh Province, Vietnam</footnote>

<footnote n="3">³ Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany</footnote>

<footnote n="4">⁴ Korea Astronomy and Space Science Institute, 776 Daedeokdae-ro, Yuseong-gu, Daejeon 34055, Republic of Korea</footnote>

<footnote n="5">⁵ University of Science and Technology, Korea, 217 Gajeong-ro, Yuseong-gu, Daejeon 34113, Republic of Korea</footnote>

<footnote n="6">⁶ Department of Astrophysics, Vietnam National Space Center, Vietnam Academy of Science and Technology, 18 Hoang Quoc Viet, Hanoi, Vietnam</footnote>

<footnote n="7">⁷ Graduate University of Science and Technology, Vietnam Academy of Science and Technology, 18 Hoang Quoc Viet, Hanoi, Vietnam</footnote>

<footnote n="8">⁸ Université de Lyon 1, ENS de Lyon, CNRS, Centre de Recherche Astrophysique de Lyon (CRAL) UMR5574, 69230 Saint-Genis-Laval, France</footnote>

<footnote n="9">⁹ School of Physics and Astronomy, University of Minnesota, 115 Union St SE, Minneapolis, MN 55455, USA</footnote>

<footnote n="10">¹⁰ Astrophysics Research Institute, Liverpool John Moores University, 146 Brownlow Hill, Liverpool L3 5RF, UK</footnote>

<footnote>* Corresponding author; nganle191919@gmail.com</footnote>
enough to prevent the gravitational collapse. The collapse of super-critical cores will lead to the formation of protostars. Conversely, when the core is sub-critical, strong B-fields can act against gravitational collapse. However, several mechanisms might enhance the mass-to-flux ratio and cause the core to collapse, including ambipolar diffusion (see e.g., Mestel & Spitzer 1956; Mouschovias et al. 2006; Liu et al. 2022) and magnetic reconnection diffusion (Lazarian & Vishniac 1999; Santos-Lima et al. 2012).

Feedback from massive stars (O- or B-type) can influence the properties of the parental molecular cloud and the star formation processes (e.g., Krumholz et al. 2014; Pabst et al. 2019, 2020). Three-dimensional magnetohydrodynamic (3D MHD) simulations have suggested that B-fields also influence the structure of clouds along with radiation pressure feedback from massive stars (e.g., Henney et al. 2009; Mackey & Lim 2011a). It has been shown that when the B-field orientations are parallel to the ionizing radiation field, the cloud will tend to have a flattened shape, while the cloud will have a column-like or pillar structure if the B-fields are perpendicular to the radiation field (Mackey & Lim 2011b). These findings from simulations have been in agreement with observations and could explain various globular and elongated substructures in the IC 1396 region (Soam et al. 2018), the famous pillar structure in the M16 region (Pattle et al. 2018), and the Horsehead Nebula (Hwang et al. 2023).

Interstellar dust grains are known to have aspherical shapes and can be systematically aligned with a fixed direction in space, for instance, interstellar B-fields (Hall 1949; Hiltner 1949). The alignment of interstellar dust grains with B-fields involves the internal alignment of the grain axis of the maximum inertia moment (also the shortest axis) with its angular momentum (i.e., internal alignment) and the alignment of the angular momentum with B-fields (external alignment). Internal alignment is caused by internal relaxation processes, including Barnett relaxation and inelastic relaxation (Purcell 1979). The external alignment is caused by the radiative torque (RAT) alignment mechanism (Draine & Weingartner 1997; Lazarian & Hoang 2007) and by the more comprehensive magnetically enhanced RAT mechanism (MRAT, Hoang & Lazarian 2016). The internal and external alignment processes establish the alignment of the grain with its minor axis parallel to the B-fields (see, e.g., Hoang et al. 2022). Dust grains absorb or scatter starlight in the ultraviolet (UV), optical, or near-infrared (NIR) range. However, aligned dust grains tend to absorb and scatter the electric field component of the starlight along their major axis more efficiently. This results in the polarization of background starlight with the polarization vector parallel to the grain's minor axis. Thus, the polarization of the starlight seen in the UV to optical or NIR wavelength has the polarization orientation parallel to the B-field orientation. Aligned dust grains also re-emit the energy absorbed from starlight in the longer wavelength from far-infrared (FIR) to sub-millimeter (sub-mm), the so-called thermal dust emission. As the electric field component of starlight is absorbed more along the grain's major axis, the dust grain will re-emit more energy with the electric field component along this major axis than along the minor axis. This causes the thermal emission of aligned dust grains to be polarized with the polarization orientation along the grain's major axis. Because dust grains have their minor axes aligned with the B-field orientation, the polarization orientation of thermal emission is perpendicular to the B-field orientation. Thus, if the position angle of the polarized thermal dust emission is known, one can rotate the polarization angle by 90° to infer the position angle of the B-fields in the medium. Linear polarized thermal dust emission has been routinely used

to trace B-fields in various environments of molecular clouds in the Milky Way (see, e.g., Ward-Thompson et al. 2017; Chuss et al. 2019; Pattle et al. 2022) and extragalactic regions (see, e.g., Lopez-Rodriguez et al. 2021, 2022; Borlaff et al. 2023). The Davis-Chandrasekhar-Fermi method (DCF; Davis 1951; Chandrasekhar & Fermi 1953) has been widely used to measure the strengths of the B-field projected on the plane of the sky, $B_{\text{pos}}$, based on the assumption that the kinetic energy of turbulent gas balances to that of turbulent B-fields.

Adopting the mean value of B-field strengths along the line of sight ($B_{\text{los}}$) measured using the Zeeman technique over different molecular clouds, Crutcher et al. (2010) derived a relationship between the B-field strengths and gas volume densities ($n_{\text{H}}$). For low densities of $n_{\text{H}} \le 300 \text{ cm}^{-3}$, B-fields were found not to considerably change with densities. However, for high densities above $n_{\text{H}} \sim 300 \text{ cm}^{-3}$, the B-field strengths increase with densities as $n_{\text{H}}^k$, where the power-law index $k$ is ~0.65. Pattle et al. (2023) compiled a large number of measurements of the B-field strength projected on the plane of the sky ($B_{\text{pos}}$) estimated using the DCF method. These authors showed that these $B_{\text{pos}}$ generally follow the results from Crutcher et al. (2010) until $n_{\text{H}} \ge 10^7 \text{ cm}^{-3}$. Beyond this $n_{\text{H}}$ value, the $B_{\text{pos}}$ distribution is more sparse (see Fig. 2a in Pattle et al. 2023). Theoretical studies have suggested that the power-law index varies between $k \sim 2/3$ for weak B-fields (Mestel 1966) and $k \le 0.5$ for strong B-fields (Mouschovias & Ciolek 1999). Therefore, the maps of the morphology and strength of the B-fields within a specific cloud are essential to understanding the role of the local B-fields on star formation within the cloud. Previous authors have attempted to derive the map of B-fields toward some well-known star-forming regions, for instance, OMC-1 (Guerra et al. 2021; Hwang et al. 2021), 30 Doradus (30 Dor, Tram et al. 2023), and the G11.11-0.12 filament (Ngoc et al. 2023). These studies revealed significant variations in the B-field morphology and strength across the clouds, demonstrating that B-fields strongly influence the cloud's evolution and, thus, star formation activities.

In this work, we focus on a specific system of the $\rho$ Ophiuchus-A (hereafter $\rho$ Oph-A) molecular cloud to study the effects of stellar feedback and B-fields on the evolution of the cloud and star formation. The $\rho$ Oph-A molecular cloud is located in the northern part of the L1688 active low-mass star-forming region in the $\rho$ Ophiuchus dark cloud complex (Loren et al. 1990; Wilking et al. 2008; Esplin & Luhman 2020). The $\rho$ Oph-A cloud was initially identified as a dense gas clump based on DCO$^+$ observations performed by Loren et al. (1990), along with other clumps named $\rho$ Oph-B to $\rho$ Oph-F. This region has high extinctions ($A_V \ge 20$ mag) corresponding to the high gas column densities compared to other regions, which are typically more quiescent (see Fig. 1). Being one of the closest low-mass star-forming regions in our Galaxy, located at a distance of ~137 pc (Ortiz-León et al. 2017), the gas and dust content in the Oph-A region has been studied extensively at almost every wavelength by both ground-based and space instruments for many years (see e.g., Wilking et al. 2008, for a review). The $\rho$ Oph-A region hosts numerous prestellar and protostellar cores identified by the Ground Belt Survey (Ladjelate et al. 2020) and protostars in the early stages (i.e., Class 0/I Young Stellar Objects (YSOs), Enoch et al. 2009; Evans et al. 2009; Connelley & Greene 2010). The $\rho$ Oph-A cloud is mostly influenced by a nearby embedded early high-mass B3/4 star (Brown & Zuckerman 1975; Mookerjea et al. 2018), Oph-S1, which has a mass of ~8 $M_\odot$ (Hamaguchi et al. 2003). The Oph-S1 star is surrounded by a compact H II region (Andre et al. 1988) and is
<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Distribution of dust temperatures (colors) and H₂ column densities (contours) in ρ Oph-A (Santos et al. 2019). Contour levels of N<sub>H</sub>₂ are shown at (1, 2, 4, 6, 8, and 10) × 10<sup>22</sup> cm<sup>-2</sup>. The star symbol shows the position of Oph-S1. Small circles show the positions of YSOs (Enoch et al. 2009; Evans et al. 2009; Connelley & Greene 2010). Yellow squares and cyan diamonds indicate the positions of prestellar and protostellar cores, respectively, identified as part of the Herschel Gould Belt Survey (Ladjelate et al. 2020). The “x” symbol shows the position of the SM1 core.</figcaption></figure>

believed to trigger the star formation activity in Oph–A, espe-
cially along the ridge of the dense cloud with indications of
protostellar clumps (see Fig. 2a in Motte et al. 1998). The star
formation activity in the ρ Oph-A is also believed to be affected
by the Sco OB2 association located in the western part of the
cloud at a distance of ~145 pc (de Zeeuw et al. 1999).

We investigated B-fields in ρ Oph-A using data observed in various wavelengths: NIR (JHK bands) using imaging polarimeter SIRPOL (Kandori et al. 2006) mounted on the 1.4-m Infrared Survey Facility (IRSF) telescope (Kwon et al. 2015), FIR by the High-resolution Airborne Wideband Camera-plus (HAWC+; Harper et al. 2018) on board the 2.7-m Stratospheric Observatory For Infrared Astronomy (SOFIA) (Santos et al. 2019), and sub-mm using the linear polarimeter POL-2 placed in front of the Submillimetre Common-User Bolometer Array 2 (SCUBA-2; Holland et al. 2013) camera mounted on the 15-m James Clerk Maxwell Telescope (JCMT) (Kwon et al. 2018). Santos et al. (2019) performed the analysis on polarimetric data from SOFIA/HAWC+ observations at 89 and 154 µm to study the change in the slope of the polarization spectra with respect to the column densities and dust temperatures in ρ Oph-A. Later, Tram et al. (2021) used a similar SOFIA/HAWC+ dataset to investigate the relation between the polarization fraction of thermal dust emission versus dust temperatures in the cloud and to understand the mechanism of dust grain alignment and disruption induced by RATs. Kwon et al. (2018) estimated the B<sub>pos</sub> values in ρ Oph-A, using the JCMT/POL-2 data at 850 µm, toward a few of the densest parts of the cloud, yielding B<sub>pos</sub> values in the range of ~0.5–5 mG. The highest strength was measured toward the densest part of the cloud associated with the starless core SM1 (see Fig. 1). These B<sub>pos</sub> strengths in ρ Oph-A are comparable to values found in the OMC–1 region (Hwang et al. 2021). Kwon et al. (2018) also found that the B-field direction is generally perpendicular to the main structure of the cloud, which is generally seen in the dense cores or high-density

structures of molecular clouds. The ρ Oph-A region provides an
ideal isolated system for studying the impact of the B-fields in
the star-forming region and investigating the interplay between
B-fields and the stellar feedback into the molecular cloud. In
this paper, we aim to: 1) determine the B-field morphology and
strengths in the entire ρ Oph-A region using the SOFIA/HAWC+
dust polarization data; 2) investigate the role of B-fields relative
to gravity and turbulence in the star formation activity in the
cloud; and 3) quantify the combined effect of B-fields and radi-
ation pressure from the nearby Oph-S1 star on the shape of the
ρ Oph-A cloud.

The paper is organized as follows. We describe a set of
archival observational data in Sect. 2. We present the results
of mapping the B-field morphology and strength using various
techniques in Sect. 3 and discuss the roles of B-fields relative
to gravity, turbulence, and radiation pressure feedback from the
Oph-S1 star in Sect. 4. Finally, we summarize our results in
Sect. 5.

## 2. Observational data

2.1. Dust polarization data from SOFIA/HAWC+

We use the FIR polarimetric observations from the onboard SOFIA/HAWC+ toward ρ Oph-A. The observations were done on 2017 May 17, under the proposal ID 70_0511 (PI: Darren Dowell), using the nod-match-chop mode, with 154° chop angle, 480" chop throw, and 4-points dithers. Here, we use the SOFIA/HAWC+ Level-4 data for two bands C and D, centered at 89 and 154 µm, respectively, which is fully calibrated and can be directly downloaded from the SOFIA archive<fnref n="1" />. The final polarimetric maps cover an area of ~4.5'×4.5' (corresponding to 0''18 × 0''18 pc²) toward the central region in ρ Oph-A with the angular resolutions of 7.8" and 13.6" full-width-half-maximum (FWHM) at 89 and 154 µm, respectively. The pixel size of HAWC+ maps is ~2.0" at 89 µm and 3.4" at 154 µm. The observations and data reduction are described in detail in Santos et al. (2019).

From the maps of the Stokes-*I*, Stokes-*Q*, and Stokes-*U* (along with their uncertainties σ_I, σ_Q, σ_U), the maps of polarization fraction (p_0) and its corresponding uncertainty (σ_p) were calculated. The de-biased polarization fraction, p, is calculated as: p^2 = p_0^2 - σ_p^2 (Vaillancourt 2006). The maps of polarization angle (θ) and its corresponding uncertainty (σ_θ) in degrees can be derived from Stokes-*Q* and Stokes-*U*. We refer to Gordon et al. (2018) for a full description of estimates of the polarization fraction, polarization angle, and their corresponding uncertainties. To enhance the quality in measuring the maps of B-field orientation and strength in ρ Oph-A, we apply some critical conditions in choosing vectors (pixels) in SOFIA/HAWC+ maps. We mask out from SOFIA/HAWC+ maps the pixels for which 1/σ_I < 54 for band C and 1/σ_I < 214 for band D. In both bands, we further removed the pixels for which we have p/σ_p < 3 and p>50% from the HAWC+ maps (see Appendix A).

2.2. Spectroscopic data

We used the spectroscopic observations of HCO⁺(4–3) line at ~356.734 GHz obtained from the Heterodyne Array Receiver Program (HARP; Buckle et al. 2009) mounted on JCMT. The observations were carried out during June and July 2011 (Proposal ID: M11AU13), covering the ρ Oph-A region with a beam

<footnote n="1">¹ https://irsa.ipac.caltech.edu/Missions/sofia.html</footnote>
<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Maps of the B-field orientation (white segments) on top of the total intensity (Stokes-*I*, color maps) toward ρ Oph-A, observed with SOFIA/HAWC+ in band C (left) and band D (right), centered at 89 and 154 µm, respectively. The gray contours in each map show the Stokes-*I* intensity at the corresponding wavelength, with contour levels at 0.12, 0.2, 0.4, 0.6, and 0.8 Jy/arcsec<fnref n="2" />. The star symbol marks the position of Oph-S1 star. Small white circle symbols indicate the positions of YSOs. The filled orange circles show the beam size of each color map.</figcaption></figure>

size (or FWHM) of ~14''. We used the reduced spectral map of HCO⁺(4–3) line directly downloaded from the JCMT science archive<fnref n="2" />. The data cube had an original spectral resolution of 0.025 km s⁻¹ and was smoothed to have a final spectral resolution of 0.5 km s⁻¹. For comparison to the HCO⁺(4–3) data, we also used observations of N₂H⁺(3–2) line at ~279.515 GHz obtained from the APEX-2 instrument on the 12-m Atacama Pathfinder EXperiment telescope (APEX; Güsten et al. 2006). The data cube has an original spectral resolution of 0.082 km s⁻¹. Here, we used the back-end data cube with a spectral resolution of 0.25 km s⁻¹, covering the ρ Oph-A region with the beam size of ~22'' (for details, see Liseau et al. 2015).

## 2.3. Dust temperature and column density maps

We adopted the maps of dust temperature (*T*<sub>dust</sub>) and gas column density (*N*<sub>H₂</sub>) from Santos et al. (2019). These maps were obtained from the fitting of a modified black body to the spectral energy distribution (SED) using the continuum flux maps from Herschel/PACS (Poglitsch et al. 2010) at 70, 100, and 160 µm, with an FWHM of 11.4'' (for details, see Sect. 3.2 in Santos et al. 2019). Figure 1 shows the spatial distribution *T*<sub>dust</sub> toward ρ Oph-A, overlaid with contours of column density, *N*<sub>H₂</sub>. In this study, we use the map of *T*<sub>dust</sub> as the map of the gas kinetic temperature, assuming that the gas is in the local thermodynamic equilibrium (LTE) and well-coupled with the dust.

# 3. Magnetic fields in ρ Oph-A

## 3.1. Magnetic field morphology

We construct the B-field morphology in ρ Oph-A by rotating polarization angles by 90°. Figure 2 shows the maps of B-field

orientation (half-vectors) inferred from maps of the polarization angle at 89 and 154 µm. The background color maps show the Stokes-*I* intensity at the respective wavelength. Here, the position angle of B-fields (θ<sub>B</sub>) is measured east of north, and in the range between 0° and 180°. The polarimetric data at both wavelengths can well reveal B-fields toward the dense region in ρ Oph-A, where Stokes-*I*>0.12 Jy/arcsec<fnref n="2" /> (the area inside the lowest contour in each map). The data at 89 µm can further resolve the B-fields morphology extended toward the warmer dust region located at the east side of the Oph-S1 star.

The B-field vectors in ρ Oph-A are generally well-ordered and largely follow the material structure in low- to intermediate-densities regimes (see Fig. 2). B-field orientations seem to be bent toward the higher density region in the center of the map and are mostly perpendicular to the ridge of the cloud. On the eastern side of the map, the orientation of B-field vectors tends to vary in the counter-clockwise direction from north to south; thus, the position angle of B-fields increases. Particularly, in the northeast part of the map, B-field vectors run along the north-east to southwest direction with a median θ<sub>B</sub> value of ~56° and ~49° at 89 and 154 µm, respectively. This feature is also seen in a small area in the northwest part of the map. Toward the east side, B-field vectors orient nearly to the east-west direction, with median θ<sub>B</sub> values of ~83° at 89 µm and ~75° at 154 µm. On the other hand, B-field vectors in the southeast region tend to orient consistently to the southeast to northwest directions, with median θ<sub>B</sub> values of ~92° and ~85° at 89 and 154 µm, respectively. Interestingly, the data at 154 µm shows B-field vectors running nearly vertically, that is, along the north-south direction toward a small area at the north part of the map, which is not seen in the data at 89 µm (see Fig. 2 and also Fig. B.1).

In addition, there are two outstanding components of B-fields toward the regions characterized by the highest H₂ column densities: one component associated with the peak of continuum emission at 154 µm runs along the northeast to southwest direction with the median B-field orientation of ~70°; the other is

<footnote n="2">² https://www.cadc-ccda.hia-iha.nrc-cnrc.gc.ca/en/jcmt/</footnote>
located at the north of the first component and associated with the peak of continuum emission at 89 µm, which exhibits the B-field vectors along the east-west direction with a median position angle of ~88° and ~81° at 89 and 154 µm, respectively (see Fig. B.1). We will further discuss the B-field morphology together with B-field strengths in Sect. 4.1.

## 3.2. Measurement of magnetic field strength using the DCF method

We first used the classical DCF method (Davis 1951; Chandrasekhar & Fermi 1953) (hereafter we refer to as the DCF method) to calculate the strengths of the B-field projected in the plane-of-the-sky ($B_{pos}^{DCF}$) in the entire region of the $\rho$ Oph-A cloud. We compared our obtained results with those obtained using the more advanced techniques described in Sect. 3.3.

According to the DCF method, the B-field strength can be estimated by Ostriker et al. (2001):

$$B_{pos}^{DCF} = f \sqrt{4\pi\rho \frac{\sigma_{NT}}{\delta\theta}}, \quad (1)$$

where $f$ is the correction factor, $\rho = \mu m_H n_{H_2}$ is the gas mass density in units of g cm$^{-3}$, $\mu=2.8$ is the mean molecular weight per hydrogen molecule (Kauffmann et al. 2008), and $m_H$ is the mass of the hydrogen atom, $n_{H_2}$ is the gas volume density in units of cm$^{-3}$. Then, $\sigma_{NT}$ is the one-dimensional (1D) non-thermal velocity dispersion along the line of sight (LOS) of the ideal gas species that traces the dense region, similarly to the dust-polarized emission in units of km s$^{-1}$, and $\delta\theta$ is the polarization angle dispersion in degrees. The correction factor, $f$, is uncertain, and numerical simulations by Ostriker et al. (2001) show $f=0.5$ when the polarization angle dispersion does not exceed 25°. However, varying values of $f$ will result in different $B_{pos}$ strengths (see e.g., Skalidis et al. 2021; Liu et al. 2021, for more detailed discussions).

Here, we use the formula given in Crutcher (2004):

$$B_{pos} \approx 9.3 \sqrt{n_{H_2}} \frac{\Delta v_{NT}}{\delta\theta} (\mu G), \quad (2)$$

where $\Delta v_{NT} = \sigma_{NT} \sqrt{8 \ln 2}$ is non-thermal FWHM velocity dispersion measured in units of km s$^{-1}$, and the correction factor $f=0.5$ has been used. This formula is widely used to estimate $B_{pos}$ via the DCF technique.

### 3.2.1. Polarization angle dispersion maps

It is essential to see how the B-field orientation changes locally within the molecular cloud. To do so, we calculated $\delta\theta$ across the $\rho$ Oph-A region. We used the method presented in Hwang et al. (2021), for which $\delta\theta$ is calculated pixel-by-pixel throughout the map of $\rho$ Oph-A. We summarize the procedure of calculating $\delta\theta$ toward the pixel $i^{th}$ in the map as follows:

1. We used a box covering the pixel $i^{th}$ with a size of 9×9 pixels, corresponding to two beam sizes of the SOFIA/HAWC+ observations (i.e., 15.6" and 36.4" in bands C and D, respectively).

2. We calculated the mean value of polarization angles within the small box, $\bar{\theta}_i$, which implies the local mean polarization angle toward the $i^{th}$ pixel.

3. For each pixel within the small box, we calculated the difference between the observed polarization angle and the mean value. We computed the root mean square (RMS) of the distribution of angle differences in the box. This RMS value

indicates the polarization angle dispersion calculated toward the $i^{th}$ pixel.

4. The uncertainty of the polarization angle dispersion toward the pixel $i^{th}$ is propagated from the uncertainty of polarization angles in the small box and determined as the mean value of polarization angle uncertainties within the small box.

By repeating the process for all pixels in the entire region, we obtained the map of $\delta\theta$ in $\rho$ Oph-A. We note that if there are more than 50% of pixels undefined (i.e., having NaN value of $\theta$) in the computing box, the result of polarization angle dispersion calculation toward pixel $i^{th}$ is ignored. We applied this criterion to keep the high quality of the $\delta\theta$ estimate. It is also necessary to de-bias the polarization angle dispersion due to the uncertainty of the observed polarized angles.

Figure 3 shows the maps of polarization angle dispersion toward $\rho$ Oph-A using HAWC+ data in two bands C (left panel) and D (right panel), respectively. The polarization angle dispersion calculated using data in bands C and D cover ranges of 2.3°–32.2° and 1.3°–31.5°, which correspond to the median values of 5.7° and 5.8°, respectively. The differences of $\delta\theta$ in two wavelengths imply that the data at the shorter wavelength tends to trace warm dust surrounding the high-mass star S1 and toward the outer layers; whereas the data at the longer wavelength probes the colder dust in the denser layers of the cloud. We note, however, that in maps of $\delta\theta$ in both bands, there are a few pixels where $\delta\theta > 25°$ (27 pixels in band C and 45 pixels in band D). Thus, we do not use the $\delta\theta$ value in these pixels to estimate $B_{pos}$ in Sect. 3.2.4.

### 3.2.2. Velocity dispersion map

We used the HCO$^+$ (4–3) molecular data observed by JCMT/HARP to estimate the total velocity dispersion, $\sigma_v$. This molecular line has a critical density of ~$10^6$ cm$^{-3}$ (Greve et al. 2009) and thus is considered a good tracer for probing the dense region in the $\rho$ Oph-A cloud. The HCO$^+$ (4–3) emission is well-represented by a single Gaussian profile. Thus, in each pixel, we fit the line with a single Gaussian. The central velocity ($v_0$) across $\rho$ Oph-A covers the range from 2.25 to 3.76 km s$^{-1}$ with a median value of 3.12 km s$^{-1}$, which is in agreement with the systemic velocity of $\rho$ Oph-A of ~3.44 km s$^{-1}$ (André et al. 2007). We used the line width of the Gaussian as the measure of total velocity dispersion toward each pixel.

To construct the map of the 1D non-thermal (turbulent) velocity dispersion ($\sigma_{NT}$), we subtracted the thermal component ($\sigma_{v,th}$) from the map of the total velocity dispersion, given by:

$$\sigma_{NT} = \sqrt{\sigma_v^2 - \sigma_{v,th}^2} = \sqrt{\sigma_v^2 - \frac{k_B T_{gas}}{m_{mol}}} (\text{km s}^{-1}), \quad (3)$$

where $m_{mol}=29$ amu is the mass of the HCO$^+$ molecule, $k_B$ is the Boltzmann constant, and $T_{gas}$ is the kinetic temperature of the gas. Here, we assume that the gas is in LTE condition and well-coupled with the dust in the cloud. Thus, the gas kinetic temperature is equal to the excitation temperature in all levels ($T_{ex}$) and equal to the dust temperature: $T_{gas} = T_{ex} = T_{dust}$. We convolved the map of $T_{dust}$ to the lower resolution of the molecular line map (14") to obtain the same pixel size on both maps. We converted the velocity dispersion of the non-thermal component into the FWHM velocity dispersion, given as $\Delta v_{NT} = \sigma_{NT} \sqrt{8 \ln 2}$. The uncertainty of non-thermal FWHMs is propagated from the uncertainty of the line width obtained from the Gaussian fit.
<figure><img src="image_5.png" /><figcaption>Fig. 3. Polarization angle dispersion maps using data in band C (left panel) and band D (right panel). The black contours are from the total intensity map at the corresponding wavelength, with the contour levels at 0.12, 0.2, 0.4, 0.6, and 0.8 Jy/arcsec². The star symbol in each map indicates the position of the Oph-S1 star.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Fig. 4. Distribution of the non-thermal FWHM of the HCO⁺(4–3) molecular line toward ρ Oph-A. The star symbol marks the position of the Oph-S1 star. White contours show the continuum emission (Stokes-I) at 154 µm, with contours levels at 0.12, 0.2, 0.4, 0.6, and 0.8 Jy/arcsec².</figcaption></figure>

Figure 4 shows the map of the 1D non-thermal FWHM of the HCO⁺ molecular line toward ρ Oph-A. The ΔvNT values are in the range from 0.62 to 2.98 km s⁻¹ with a median value of 1.19 km s⁻¹. We find a good agreement between our result obtained with HCO⁺(4–3) and those from other line tracers, namely, 0.95 km s⁻¹ using N₂H⁺(3–2) from APEX and 1.5 km s⁻¹ using C¹⁸O(3–2) from JCMT/HARP (see Appendix C for details). Thus, it is plausible to use results from the HCO⁺(4–3) line in the following analysis.

### 3.2.3. Gas volume density map

We estimated the gas volume, $n_{H_2}$, across the dense cloud from $N_{H_2}$. Here, we assume that the cloud has a spherical

Plummer-like shape (Liseau et al. 2015; Santos et al. 2019), which is frequently used to describe the globular dense cores (Plummer 1911). We consider a sphere centered at the peak position of $N_{H_2}$, located at (RA, Dec)= (16h26m27s8, -24°24'00"0), which is also associated with the position of the starless core SM1 (the separation is much less than 1"). We took the radius of the cloud $R_{cl} \approx 120''$ (~0.074 pc, Santos et al. 2019). To describe the effect of stellar feedback from the high-mass Oph-S1 star, we further add a spherical shell with a radius of $R_{S1} \approx 80''$ surrounding this star. We also assume that the positions of peak column density and the S1 star are located in the same plane perpendicular to the LOS. The top panels in Fig. 5 show the 3D schematic view of the cloud (left) and a cross-section on the xz-plane at a given y position (right). The x- and y-axes perpendicular to each other are in the plane of the sky and the z-axis is along the LOS. The center of the sphere is located at $(x_0, y_0, z_0) = (0, 0, 0)$. The depth of the cloud at each position (x, y) along z-axis is calculated by:

$$W(x, y) = 2(z_1 - z_2) = 2\left(\sqrt{R_{y}^{2} - x^{2}} - \sqrt{R_{y,S1}^{2} - (x - x_{S1})^{2}}\right), \quad (4)$$

where $R_y = \sqrt{R_{cl}^2 - y^2}$ and $R_{y,S1} = \sqrt{R_{S1}^2 - (y - y_{S1})^2}$ are the radius of horizontal circles perpendicular to the plane of sky at y position (see top right panel in Fig. 5). We obtain the map of W, shown in the bottom left panel in Fig. 5. The values of W are in the range of ~20''–222''. We then infer $n_{H_2}$ from $N_{H_2}$ using the formula:

$$n_{H_2} = \frac{N_{H_2}}{W}. \quad (5)$$

The bottom right panel in Fig. 5 shows the spatial extent of $n_{H_2}$ across ρ Oph-A, overlaid by the contours representing the continuum emission in band D (Stokes-I). The $n_{H_2}$ map has an angular resolution of 11.4" and a pixel size of ~3.2",
<figure><img src="image_7.png" /><figcaption>**Fig. 5.** Cartoon of 3D structure of the ρ Oph-A cloud where we assume that the cloud has a spherical shape but is pushed by radiation feedback from a nearby high-mass star Oph-S1, shown at the top (left panel, see text) and a plane of the constant coordinate y with the LOS $\hat{z}$ (right panel). Bottom row displays the map of the depth of the cloud, *W* (left panel), and the map of the gas volume density ($n_{H_2}$) inferred from $N_{H_2}$ (right panel). In each panel, the “×” and star symbols mark the peak position of $N_{H_2}$ and location of Oph-S1 star, respectively. Toward maps of *W* and $n_{H_2}$ in bottom panels, the black dashed circle has a radius $R_{cl}=120''$ and black solid contours show the continuum emission at 154 µm with levels similar as shown in Fig. 4.</figcaption></figure>

similar to that of the $N_{H_2}$ map. The gas density spans from ~$10^{4.4} – 10^6$ cm$^{-3}$ and is highest at the center of the cloud. The pattern of $n_{H_2}$ distribution well follows the intensity map in band D, revealing the ridge shape of the dense cloud. Some pixels located at the circular boundary in the $n_{H_2}$ map show anomalously high densities; this is due to the very low depths caused by the spherical model, and thus we exclude these pixels from our $B_{pos}$ strength measurement.

### 3.2.4. Magnetic field strength maps

We used the maps of $\delta\theta$, $\Delta v_{NT}$, and $n_{H_2}$ to measure $B_{pos}^{DCF}$ in the entire region of ρ Oph-A. All input maps were smoothed to the angular resolution of 14", corresponding to that of the $\Delta v_{NT}$ map. We then reconstructed the maps of $B_{pos}$ at 89 and 154 µm toward ρ Oph-A using Eq. (2), which are shown in Fig. 6.

The $B_{pos}^{DCF}$ values vary from 0.2–2.3 mG and from 0.2–2.5 mG for data at 89 and 154 µm, respectively. The highest

strength is found at the center region of the cloud, associated with the emission peak at 154 µm and also the position of the SM1 core. The median $B_{pos}^{DCF}$ values are 0.7 ± 0.3 and 0.7 ± 0.4 mG<fnref n="3" />, as measured using 89 and 154 µm data, respectively. The uncertainty of $B_{pos}$ is propagated based on the uncertainty of the velocity dispersion of the gas and polarization angle dispersion. Since $n_{H_2}$ map is inferred from $N_{H_2}$ resulted in SED fitting models (see Santos et al. 2019). In addition, the estimated depth of the cloud is based on the 3D geometrical shape of the cloud and also is uncertain. Thus, we did not take into account the uncertainty of the gas volume density in the total uncertainty of the B-field strengths. For this reason, the uncertainty of $B_{pos}$ might be underestimated.

<footnote n="3">³ The uncertainties of these median values are solely of the estimated $B_{pos}^{DCF}$ distribution across the cloud.</footnote>
<figure><img src="image_8.png" /><figcaption>Fig. 6. Maps of B-field strengths projected in the plane-of-the-sky ($B_{pos}$) in $\rho$ Oph-A using data in band C (left panel) and band D (right panel), respectively. The contours in each map are similar as in Fig. 2. The star symbol marks the position of the Oph-S1 star and small circle symbols show the position of YSOs.</figcaption></figure>

### 3.3. Measurement of magnetic field strength using other methods

Several new techniques have been proposed to improve the DCF method for the last decades. For instance, Hildebrand et al. (2009) introduced the angle dispersion function (ADF) method where the dispersion in polarization angles $\delta\theta$ in Eq. (2) is replaced by the square root of the structure function of polarization angles. This provides an improved measurement of B-field strength but also gives us the ratio of the turbulent to large-scale B-field strengths. Alternatively, Skalidis & Tassis (2021) introduced a modification of the DCF method (hereafter ST), assuming that the mean (ordered) and fluctuation components of the total B-field energy, $B_0$ and $B_t$, respectively, are not independent. More recently, the differential measure approach (DMA; Lazarian et al. 2022) is considered promising to enhance the accuracy of B-field strength measurements at scales of $\ell$ smaller than the turbulence injection scale and takes into account the anisotropic properties of MHD turbulence (given by factor $f'$). For the sake of cross-checking, the B-field strength measurement in $\rho$ Oph-A using the DCF method, we re-calculated the mean B-field strengths using the ADF, ST, and DMA methods at the two subregions (a) and (e) associated with the densest parts of the cloud (these subregions are named as in Kwon et al. (2018), see also in Fig. B.1).

#### 3.3.1. ADF method

According to the ADF method, the $B_{pos}^{ADF}$ can be estimated using Eq. (2) but the polarization angle dispersion $\delta\theta$ is replaced by the structure function of the polarization angles, so-called angular structure function (Hildebrand et al. 2009). The angular structure function, $D_{\theta}^{1/2}$, at an angular scale $\ell$ is estimated as:

$$D_{\theta}^{1/2}(\ell) \equiv \langle \Delta\theta^2(\ell) \rangle^{1/2} = \left\{ \frac{1}{N(\ell)} \sum_{i=1}^{N(\ell)} \left[ \theta(x) - \theta(x + \ell) \right]^2 \right\}^{1/2}, \quad (6)$$

where $\langle...\rangle$ denotes an average, $\Delta\theta(\ell) \equiv \theta(x) - \theta(x + \ell)$ is the polarization angle difference between individual pairs of polarization vectors at position $x$ and $x + \ell$, where $\theta(x)$ and $\theta(x + \ell)$ are the corresponding polarization angles at position $x$ and $x + \ell$, respectively, and $N(\ell)$ is number of pairs of vectors with a displacement of $\ell$. At scales $\ell$ much smaller than the scale for a variation of the large-scale B-field structure (typically of ~5' seen by Planck observations), the total angular dispersion function can be expressed as:

$$D_{\theta}(\ell) \approx b^2 + m^2 \ell^2 + \sigma_{M}^2(\ell), \quad (7)$$

where $b$ represents the turbulent component of B-fields ($B_t$), $m\ell$ describes the large-scale structure of B-fields ($B$), and $\sigma_M(\ell)$ is subject to measurement uncertainties. Also, the ratio between the turbulent and large-scale B-field strength indicates the angular dispersion, given as (Hildebrand et al. 2009):

$$\delta\theta \approx \frac{\langle B_t^2 \rangle^{1/2}}{B} = \frac{b}{\sqrt{2-b^2}}. \quad (8)$$

When $B_t \ll B$, i.e., $b \ll 1$ rad, we can obtain $\delta\theta \sim b/\sqrt{2}$.

For each subregion (a) and (e) in $\rho$ Oph-A, we calculate the angular structure function $D_{\theta}(\ell)$ using HAWC+ data at 154 µm via Eq. (6) and then perform the fit to the quadratic function described in Eq. (7) using the Python curve_fit routine. We select data points that $13.6'' < \ell \le 50''$ to avoid fitting to structures higher than the beam size and to better fit to smaller-scale structures. Top panels in Fig. 7 show structure functions of polarization angles as a function of the distance of displacement $\ell$ toward subregions (a) and (e). From the fitted $b$ value of the fits, we derive $\delta\theta$ values of $10.9^\circ \pm 0.2^\circ$ and $17.5^\circ \pm 0.2^\circ$ toward subregions (a) and (e), respectively. We find these values in good agreement with mean polarization angle dispersion values estimated using our approach described in Sect. 3.2.1, which are ~9.2° and ~14.3° in subregions (a) and (e), respectively. Toward subregions (a) and (e), mean gas volume densities
<figure><img src="image_9.png" /><figcaption>Fig. 7. Structure function of polarization angles (top panels) and of central velocity (middle panels) measured in the subregions (a) (left) and (e) (right). The bottom panels show the resulting $B_{pos}^{DMA}$ strength estimated via Eq. (10) in the two subregions (a) and (e). The HAWC+ data at 154 µm is used to calculate the structure-function (data points). The red solid line in the top panels shows the fitted line of the angular structure-function to the quadratic function, where we chose data points with $13.6'' < l \le 50''. The blue dash-dotted line in each panel indicates the upper $l$-limit of $50''$ used for fitting the angular structure-function and for estimating the mean value of $B_{pos}^{DMA}$ in the two subregions (a) and (e) (see texts). The blue dashed line in each panel indicates the beam size of HAWC+ observation at 154 µm.</figcaption></figure>

are $2.1 \times 10^5$ and $1.1 \times 10^5$ cm$^{-3}$ respectively, and their mean non-thermal FWHMs are 1.74 and 1.45 km s$^{-1}$, respectively. We substitute new $\delta\theta$ values together with the above mean values of $n_{H_2}$ and $\Delta\nu_{NT}$ into Eq. (2) to calculate $B_{pos}^{ADF}$ strengths in subregions (a) and (e), which are 0.69 mG and 0.26 mG, respectively. These values are in good agreement with $B_{pos}^{DCF}$ strengths obtained using the DCF method, where $B_{pos}^{DCF}$ values are in ranges of 0.32–1.55 mG and 0.15–0.97 mG toward subregions (a) and (e), respectively, with corresponding mean $B_{pos}^{DCF}$ values of ~0.84 and ~0.37 mG.

### 3.3.2. ST method

The ST method is applicable for the compressible B-fields mode and thus is expected to provide a more accurate estimate of B-field strengths (Skalidis et al. 2021). Another advantage of the ST method is that this method accounts for the correlation between the fluctuation and mean components of the total B-field energy not to be zeros. Here, we measure maps of $B_{pos}$ strength ($B_{pos}^{ST}$) in band C and band D using the formula as appeared in

Skalidis & Tassis (2021):

$$B_{pos}^{ST} = \sqrt{4\pi\rho\frac{\Delta\nu_{NT}}{\sqrt{2\delta\theta}}} = \frac{1}{\sqrt{2\delta\theta}}B_{pos}^{DCF} \quad (9)$$

We show the maps of $B_{pos}^{ST}$ across $\rho$ Oph-A in band C and band D in Fig. B.3. We find that $B_{pos}^{ST}$ value ranges from 0.03–0.7 mG at 89 µm and 0.02–1.3 mG at 154 µm, appeared to be lower than $B_{pos}^{DCF}$. A median $B_{pos}^{ST}$ value of ~0.2 ± 0.1 mG is seen in both bands. For a direct comparison with those using other methods, we estimate the mean value of $B_{pos}^{ST}$ at 154 µm toward two subregions (a) and (e), which are ~0.19 and ~0.07 mG, respectively. These values are about 4.5 times lower than those estimated using the DCF method.

### 3.3.3. DMA method

We also calculate $B_{pos}$ strengths in two subregions (a) and (e) using the DMA method, $B_{pos}^{DMA}$, using Eq. (10) following
Lazarian et al. (2022):

$$B_{\text{pos}} = f' \sqrt{4\pi\rho} \frac{D_{\theta}^{1/2}(\ell)}{D_{\theta}^{1/2}(\ell)}, \qquad (10)$$

where $D_{\theta}^{1/2}(\ell)$ and $D_{\theta}^{1/2}(\ell)$ are the second order structure function of polarization angles and central velocity at an angular scale $\ell$, $f'$ is the correction factor depending on the anisotropic properties of MHD turbulence. Here, we adopt $f'=1$ related to the naive formalism as suggested in Lazarian et al. (2022).

We calculate the structure function of the central velocity of HCO$^{+}$(4–3) line toward subregions (a) and (e), using the similar form in Eq. (6) but for the central velocity $v_0$, shown in middle panels in Fig. 7. We then measure $B_{\text{pos}}^{\text{DMA}}$ strengths in two subregions (a) and (e) via Eq. (10), using the information from the $D_{\theta}^{1/2}(\ell)$ together with $D_{\theta}^{1/2}(\ell)$ and the mean values of $n_{\text{H}_2}$ in the subregions. The $B_{\text{pos}}^{\text{DMA}}$ as a function of displacement $\ell$ is shown in bottom panels of Fig. 7 for subregions (a) and (e). It can be seen that the resulting $B_{\text{pos}}^{\text{DMA}}$ in the range $13.6'' < \ell < 50''$ is likely flatten within their uncertainty. Thus, in each subregion we take the mean of $B_{\text{pos}}^{\text{DMA}}$ in this range to be the average value, which are ~0.34 mG and ~0.12 mG in subregions (a) and (e), respectively. These values are slightly higher than the resulting $B_{\text{pos}}$ estimated using the ST method, but they are at the low-end of that measured using the DCF method (i.e., ~2–3 times less than the mean values estimated using the DCF method presented in Sect. 3.2). These differences are indeed expected since the DCF method somewhat overestimates the $B_{\text{pos}}$ strength in the cloud. We discuss this issue further in Sect. 4.1.

### 3.4. Mass-to-flux ratio

The relative importance of B-fields compared to gravity in the cloud can be quantified by the mass-to-flux ratio, $M/\phi$, which can be estimated in units of the critical ratio, $(M/\phi_B)_{\text{crit}} = 1/(2\pi\sqrt{G})$ (Nakano & Nakamura 1978): $\lambda \equiv (M/\phi)_{\text{observed}} / (M/\phi)_{\text{crit}}$. We measure $\lambda$ in the entire region of $\rho$ Oph-A, using the formula introduced by Crutcher (2004):

$$\lambda = 7.6 \times 10^{-21} \frac{N_{\text{H}_2}}{B}, \qquad (11)$$

where $N_{\text{H}_2}$ and $B$ are in units of cm$^{-2}$ and µG, respectively.

Due to the lack of $B_{\text{los}}$ component measurements in the $\rho$ Oph-A cloud, we adopted the $B_{\text{pos}}$ component as the total B-field strengths $B$: $B \equiv B_{\text{tot}} \approx B_{\text{pos}}$. For the analysis, we employed only the $B_{\text{pos}}^{\text{DCF}}$ map estimated using the DCF method and inferred from the polarimetric data in band D; this is because it covers a larger area of the higher density region in $\rho$ Oph-A than that coming from the data in band C.

The top panel in Fig. 8 shows the spatial distribution of the mass-to-flux ratio in $\rho$ Oph-A. It reveals an increasing gradient of $\lambda$ from the central and eastern parts of the clump to the outer regions, spanning from 0.03–2.56, with a median value of ~0.26. The region exhibiting the lowest values of $\lambda \ll 1$ is close to the Oph-S1 star, suggesting that the region is firmly magnetically sub-critical, namely, B-fields are strong enough to prevent the gravitational collapse. This can also be explained by the influence of the strong winds and radiation fields from the Oph-S1 star which tend to push the material away to its vicinity at the west side. The mass-to-flux ratio increases toward the outskirt regions of the Oph-S1 star, especially toward regions along the ridge of the dense cloud (i.e. center, north, and southeast parts

<figure><img src="image_10.png" /><figcaption>Fig. 8. Maps of $\lambda$ (top panel), Alfvén Mach number (middle panel), and plasma $\beta$ parameter (bottom panel) in $\rho$ Oph-A using the $B_{\text{pos}}^{\text{DCF}}$ in band D. In each map, yellow square and cyan diamond symbols indicate the positions of prestellar and protostellar cores, respectively. Green circles indicate the positions of YSOs (Enoch et al. 2009; Evans et al. 2009; Connelley & Greene 2010). The star and “x” symbols mark the positions of the Oph-S1 star and the SM1 core, respectively. Black contours show $N_{\text{H}_2}$, similar as in Fig. 1.</figcaption></figure>

of the map) where the gas column density is $\ge 4 \times 10^{22}$ cm$^{-2}$. These regions are magnetically super-critical and so, the self-gravitational collapse might have occurred there, resulting in several protostars and a single protostellar core identified in the region (Enoch et al. 2009; Evans et al. 2009; Connelley & Greene 2010).
### 3.5. Alfvénic Mach number

To investigate the interplay between B-field and turbulence due to the non-thermal motion of the gas in $\rho$ Oph-A, we calculate the map of Alfvénic Mach number, given as:

$$ \mathcal{M}_A = \sqrt{3} \frac{\sigma_{NT}}{v_A}, \qquad (12) $$

where $\sigma_{NT}$ is given by Eq. (3), and $v_A$ is the Alfvénic speed driven by the B-fields ($v_A = B/\sqrt{4\pi\rho}$), both in units of km s$^{-1}$. The factor $\sqrt{3}$ in Eq. (12) assumes that turbulence is isotropic, and allows the calculation of the 3D Alfvénic Mach numbers. We note that adding this factor will provide upper limits to the estimated $\mathcal{M}_A$.

The middle panel in Fig. 8 shows the map of $\mathcal{M}_A$ toward $\rho$ Oph-A using results of polarimetric data at 154 µm and HCO$^+$ (4–3) as the gas tracer. The $\mathcal{M}_A$ values are in the range from 0.1–1.5 across the cloud, with a median value of 0.38. Most of $\rho$ Oph-A shows $\mathcal{M}_A < 1$ (magnetically sub-Alfvénic), suggesting the B-field dominates turbulence in the cloud. However, two small regions at the north and southwest of the map exhibit $\mathcal{M}_A \ge 1$, indicating strong Alfvénic motions of the gas (magnetically super-Alfvénic). These regions are also associated with the positions of a few prestellar and protostellar cores, and YSO candidates (Enoch et al. 2009; Evans et al. 2009; Connelley & Greene 2010).

### 3.6. Plasma $\beta$ parameter

To compare the magnetic pressure ($P_{mag}$) and thermal pressure ($P_{th}$) in $\rho$ Oph-A, we constructed the map of the plasma beta parameter ($\beta$), given as:

$$ \beta = \frac{P_{th}}{P_{mag}} = 2 \left( \frac{c_s}{v_A} \right)^2, \qquad (13) $$

where $c_s = \sqrt{k_B T_{gas}/(\mu m_H)}$ is the thermal sound speed measured in units of km s$^{-1}$, assuming $T_{gas} = T_{dust}$ (see Sect. 3.2.2).

The bottom panel in Fig. 8 shows the map of the $\beta$ parameter toward $\rho$ Oph-A in band D. We obtained $\beta < 1$ across the map with a median value of 0.02. This result implies that the cloud is supported by strong magnetic pressure, which is much more prominent than that induced by the thermal motions of the gas.

### 3.7. Virial mass analysis

To probe the dynamical state of the $\rho$ Oph-A cloud, we calculated the virial parameter ($\alpha_{vir}$), which is the ratio of the virial mass ($M_{vir}$) and the isothermal mass ($M_{iso}$) of the cloud. Here, we assume the cloud as a uniform density sphere with an effective radius of $R_{eff}$ as a whole, having a mean mass density of $\langle\rho\rangle$. The isothermal mass, $M_{iso}$, considering the self-gravity of the cloud can be estimated as:

$$ M_{iso} = \frac{4}{3} \pi \langle\rho\rangle R_{eff}^3. \qquad (14) $$

The virial mass, $M_{vir}$, determined from magnetic, thermal, and non-thermal kinetic energies can be estimated as (see e.g., Bertoldi & McKee 1992; Pillai et al. 2011):

$$ M_{vir} = \frac{5R_{eff}}{G} \left( \frac{\langle v_A \rangle^2}{6} + \langle c_s \rangle^2 + \langle \sigma_{NT} \rangle^2 \right), \qquad (15) $$

where $\langle v_A \rangle$, $\langle c_s \rangle$, and $\langle \sigma_{NT} \rangle$ are the mean values of the Alfvénic velocity, thermal sound speed, and of non-thermal velocity dispersion toward the cloud, respectively, while $G$ is the gravitational constant.

We calculate $R_{eff}$ by selecting an ellipse covering the area where $B_{pos}$ was measured. This ellipse is centered at (RA, Dec) = (16h26m29s5, -24°23'50.0') with the semi-major and semi-minor axes are $a=120''$ and $b=75''$, respectively. The $R_{eff}$ is then determined as $R_{eff} = \sqrt{ab} \approx 95''$ (~0.063 pc assuming the distance to $\rho$ Oph-A of 137 pc). Within this area, the mean gas volume density is $\langle n_{H_2} \rangle = 2.2 \times 10^5$ cm$^{-3}$, corresponding to the mean mass density $\langle\rho\rangle = 1.01 \times 10^{-18}$ g cm$^{-3}$. The corresponding mean values of velocities are: $\langle v_A \rangle = 3.36$ km s$^{-1}$, $\langle c_s \rangle = 0.09$ km s$^{-1}$, and $\langle \sigma_{NT} \rangle = 0.59$ km s$^{-1}$. We used these mean values of the velocities, non-thermal velocity dispersion, and the gas volume density to calculate $M_{iso}$ and $M_{vir}$ in Eqs. (14) and (15), yielding $\alpha_{vir} = M_{vir}/M_{iso} \approx 15$. This shows that the cloud is gravitationally unbound and that B-fields and turbulence can support the cloud against self-gravity. This high value of $\alpha$ obtained in the cloud overalls could explain the low efficiency of star formation in $\rho$ Oph-A; however, in a few regions, the density is very high (up to $10^5$ or $10^6$ cm$^{-3}$). This is consistent with the lack of prestellar cores and protostellar sources in the central region, which can be seen in the lower density parts of the cloud only (see e.g., Fig. 8). To study smaller scales across the clump, polarimetric observations with a higher spatial resolution would be required.

## 4. Discussion

### 4.1. Morphology and strengths of magnetic field in $\rho$ Oph-A

The B-field orientation in $\rho$ Oph-A is well-ordered and predominantly perpendicular to the ridge of the dense cloud (Sect. 3.1). The feature is illustrated most clearly by the data at the longer wavelength (154 µm). The data at 89 µm better traces warmer dust in the lower density region of $\rho$ Oph-A, which surrounds the Oph-S1 star.

The B-field properties measured using SOFIA/HAWC+ polarization data are qualitatively consistent with previous measurements of B-field orientations in $\rho$ Oph-A obtained at 850 µm using SCUBA-2/POL-2 (Kwon et al. 2018) and SCUBA/SCUPOL (Matthews et al. 2009) on JCMT and at NIR wavelengths (JHK bands) using IRSF/SIRPOL (Kwon et al. 2015). In the POL-2 dataset, $\rho$ Oph-A was divided into ten subregions, identified based on their B-field position angle and the polarization fraction, covering the entire extent of HAWC+ maps (see Fig. B.1). The highest density region is characterized by an average B-field position angle of 50°, which is also consistent with the NIR data toward the L1688 main cloud (Kwon et al. 2015). Table 1 shows the median B-field position angle ($\theta_B$) in each subregion, which has been covered by the FIR and sub-mm observations a well. The median B-field position angles from the 89, 154, and 850 µm datasets are in rather good agreement, with the difference in $\theta_B$ values ranging from 3° to 19° (≤30%). These differences could be explained by the fact that the sub-mm data probe colder dust layers of the cloud, while the FIR data trace warmer layers. Indeed, such foreground layers were identified through the detection of the self-absorption features in the spectra of the [O I] line at 63 µm (Liseau et al. 2006) and the [C II] line at 158 µm (Mookerjea et al. 2018).

Kwon et al. (2018) calculated $B_{pos}$ in three subregions associated with the densest parts of $\rho$ Oph-A using non-thermal velocity dispersion from the $N_2H^+$ (1–0) line (André et al. 2007)
<table><caption>**Table 1.** B-field position angle ($\theta_B$) and strength ($B_{pos}$) in subregions for the dataset at 89, 154, and 850 µm.</caption><thead><tr><th>Subregion<sup>(a)</sup></th><th>θ<sub>B,89</sub></th><th>θ<sub>B,154</sub><br>(deg)</th><th>θ<sub>B,850</sub></th><th>θ<sub>B,89–850</sub></th><th>Diff.<sup>(b)</sup><br>(%)</th><th>θ<sub>B,154–850</sub></th><th>Diff.<sup>(b)</sup><br>(%)</th><th>B<sub>pos,89</sub></th><th>B<sub>pos,154</sub><br>(mG)</th><th>B<sub>pos,850</sub></th></tr></thead><tbody><tr><td>(a)</td><td>73.9 ± 3.1</td><td>70.8 ± 3.2</td><td>54.4 ± 1.5</td><td>36</td><td></td><td>30</td><td></td><td>4.5</td><td>2.1</td><td>5.0</td></tr><tr><td>(b)</td><td>56.3 ± 2.7</td><td>48.7 ± 0.6</td><td>39.7 ± 4.4</td><td>42</td><td></td><td>23</td><td></td><td>0.1</td><td>0.2</td><td>—</td></tr><tr><td>(e)</td><td>87.9 ± 2.4</td><td>80.6 ± 3.4</td><td>99.6 ± 2.7</td><td>12</td><td></td><td>19</td><td></td><td>0.7</td><td>0.3</td><td>0.8</td></tr><tr><td>(f)</td><td>60.1 ± 1.5</td><td>55.0 ± 1.2</td><td>76.2 ± 4.7</td><td>21</td><td></td><td>28</td><td></td><td>0.2</td><td>0.3</td><td>—</td></tr><tr><td>(g)</td><td>82.9 ± 4.2</td><td>75.4 ± 1.4</td><td>66.3 ± 3.3</td><td>25</td><td></td><td>14</td><td></td><td>0.5</td><td>0.7</td><td>—</td></tr><tr><td>(i)</td><td>92.0 ± 4.3</td><td>84.9 ± 1.6</td><td>75.0 ± 3.0</td><td>23</td><td></td><td>13</td><td></td><td>0.4</td><td>0.5</td><td>—</td></tr></tbody></table>

**Notes.** (a) The subregions are named as in Kwon et al. (2018). (b) The relative difference in $\theta_{B,\lambda-850}$ is taken as $|1 - \theta_{B,\lambda}/\theta_{B,850}| \times 100$, where $\theta_{B,\lambda}$ is the B-field position angle at wavelength $\lambda = [89, 154]$ (µm).

<table><caption>**Table 2.** Comparison of $B_{pos}$ strength, mass-to-flux ratio ($\lambda$), Alfvénic Mach number ($M_A$), and $\beta$ plasma toward subregions (a) and (e) measured using different methods.</caption><thead><tr><th rowspan="2">Subregion</th><th colspan="4">B<sub>pos</sub> (mG)</th><th colspan="4">λ</th><th colspan="4">M<sub>A</sub></th><th colspan="4">β</th></tr><tr><th>DCF</th><th>ADF</th><th>ST</th><th>DMA</th><th>DCF</th><th>ADF</th><th>ST</th><th>DMA</th><th>DCF</th><th>ADF</th><th>ST</th><th>DMA</th><th>DCF</th><th>ADF</th><th>ST</th><th>DMA</th></tr></thead><tbody><tr><td>a</td><td>0.84</td><td>0.69</td><td>0.19</td><td>0.34</td><td>0.91</td><td>1.12</td><td>4.03</td><td>2.49</td><td>1.57</td><td>1.92</td><td>6.93</td><td>3.88</td><td>0.02</td><td>0.03</td><td>0.37</td><td>0.12</td></tr><tr><td>e</td><td>0.37</td><td>0.26</td><td>0.07</td><td>0.12</td><td>0.84</td><td>1.22</td><td>4.20</td><td>1.99</td><td>0.62</td><td>0.90</td><td>3.09</td><td>1.96</td><td>&lt;0.01</td><td>&lt;0.01</td><td>0.11</td><td>0.04</td></tr></tbody></table>

and the sub-mm data at 850 µm. In our study, we calculated the pixel-by-pixel $B_{pos}$ toward the entire $\rho$ Oph-A region. To facilitate comparisons with results obtained from the sub-mm data, we calculated the mean value of $B_{pos}$ in two out of these three subregions, covered by the HAWC+ data (see Table 1). Mean values of $B_{pos}$ using FIR data are lower than that obtained with the sub-mm data, which traces colder dust layers in the cloud characterized also by the higher densities. The HAWC+ data, on the other hand, essentially traces the warmer layers in the lower density region. In addition, Kwon et al. (2018) adopted values of the polarization angle dispersion as the standard deviation of the polarization angles, which are much smaller than those using the spread of the dust polarization angles (see Sect. 3.2.1). As a consequence, the $B_{pos}$ obtained with sub-mm data is expected to be stronger than the one obtained with FIR data. We find that our $B_{pos}$ estimate for $\rho$ Oph-A is typically higher than those in other molecular clouds, ranging from 50 to 400 µG (for e.g., Coudé et al. 2019; Ngoc et al. 2021; Pattle et al. 2021; Ward-Thompson et al. 2023), but similar to those found in the DR21 (OH) region (Poidevin et al. 2013) or the OMC-1 cloud in Orion A region (Pattle et al. 2017; Hwang et al. 2021) where $B_{pos}$ of a few mG was measured.

In this study, we adopt $B_{pos}$ component as the total B-field strength. If we consider the contribution of the LOS component, the total B-field strength will be stronger than what we measured. The mass-to-flux ratio, Mach Alfvénic number, and plasma $\beta$ parameter will decrease, but this will not significantly change our conclusion that B-fields are generally dominant over gravity, turbulent, and gas kinetic energies, especially toward the center of the cloud. Furthermore, as presented in Sect. 3.3, we test that the $B_{pos}$ strengths obtained with the DCF method in subregions (a) and (e) are rather in good agreement with values obtained with the ADF method but higher than those using the ST and DMA methods (see Table 2). We use new $B_{pos}$ values estimated with ADF, ST, and DMA methods toward two subregions (a) and (e) to estimate the mass-to-flux ratio, Mach Alfvénic number, and $\beta$ plasma, which are all higher than that measured with the DCF method (see Table 2). Subregion (a)

changes from magnetically sub/trans-critical to super-critical but remains to be super-Alfvénic, and has $\beta < 1$ in overall. A similar behavior is seen in subregion (e) except that this subregion changes from trans- to super-Alfvénic (see Table 2).

We go on to assume that the B-field strength in the entire cloud decreases by a factor of ~4.5, as measured with the ST method. In this case, the majority of the cloud will become magnetically trans- or super-critical ($\lambda \ge 1$), except for the regions in the vicinity of the Oph-S1 star (see the top panel of Fig. D.2). However, this trend of $\lambda$ remains similar as obtained using $B_{pos}^{DCF}$. The cloud will become mostly super-Alfvénic with $M_A > 1$, except for the southeast region of the cloud (see the middle panel of Fig. D.2). This assumption will also result in the Alfvénic velocity decreasing by a factor of ~4.5 and thus, $M_{vir} \approx 2M_{iso}$. This implies that the cloud is in the magneto hydrostatic equilibrium and may undergo gravitational collapse due to turbulent fluctuations on smaller scales.

## 4.2. Magnetic field orientation versus density structure of $\rho$ Oph-A

To investigate what is the relative orientation of the B-fields with respect to the gas structure in $\rho$ Oph-A, we quantify the projected offset angle between the orientations of B-field vectors projected in the plane-of-the-sky and the main structure associated with the ridge of the cloud, given as $\Delta\theta = \theta_B - \theta_{ridge}$, where $\theta_{ridge}$ is the position angle of the ridge. We identified the shape of the ridge using the Python package RadFil<fnref n="4" /> (Zucker & Chen 2018, see Appendix E).

Figure 9 shows the map of offset angles ($\Delta\theta$) in $\rho$ Oph-A, overlaid by contours of continuum emission at 154 µm. Clearly, the offset angle appears to increase when moving from the lower density regions to the denser regions, indicating that the relative orientation of B-fields with respect to the main structure of the cloud turns from parallel to perpendicular toward the high-density region. This is in agreement with Planck polarization

<footnote n="4">⁴ https://github.com/catherinezucker/radfil</footnote>
<figure><img src="image_11.png" /></figure>

**Fig. 9.** Map of the offset angle between B-field orientations and main cloud structure in ρ Oph-A. White contours show the continuum emission at 154 µm with similar levels as in Fig. 4. The dashed circles centered at the position of the Oph-S1 star have a radius of 60'' (black), 90'' (blue), and 120'' (green), corresponding to the levels of $A_V$ of ~6, 30, and 200 mag. The star symbol marks the position of the Oph-S1 star.

observations at 353 GHz (850 µm) toward the entire L1688 region where the relative orientation of B-fields, with respect to the gas structure changes gradually. It runs parallel in the low-density region to perpendicular toward the higher density region, with the transition occurring at $N_{H_2} \sim 3 \times 10^{21} \text{ cm}^{-2}$ and corresponding to the visual extinction $A_V$ of ~3.3 mag<fnref n="5" /> (Soler 2019). This transition has also been observed in other nearby molecular clouds in the Gould Belt using Planck 353 GHz polarization observations, typically occurring at $N_{H_2} \sim 2.5 \times 10^{21} \text{ cm}^{-2}$ or $A_V \sim 2.7$ mag (see e.g., Planck Collaboration Int. XXXV 2016; Soler et al. 2017). This phenomenon can be explained by the motion of the high-density material during the self-gravitational collapse. When the contraction occurs, the material becomes denser and compressed, which could cause the B-field lines to be re-configured. Eventually, the B-field lines are dragged and reshaped to become more uniform and perpendicular to the high-density structure of the field.

The HAWC+ data traces B-fields in the dense region of ρ Oph-A with $N_{H_2} \sim (3–240) \times 10^{21} \text{ cm}^{-2}$ (see Fig. 1), which corresponds to $A_V$ in range of ~3.2–256 mag. Thus, we further investigate in detail the relation between the relative orientation of B-fields with respect to the cloud structure and the $A_V$ within ρ Oph-A, as shown in Fig. 10. Clearly, in the regime where $A_V \le 30$ mag the offset angle ($\Delta\theta$) decreases with $A_V$ and drops from ~60° at $A_V=4$ mag down to 0° at $A_V \sim 20$ mag. This indicates that the relative orientation of B-fields changes from perpendicular back to parallel at $A_V \sim 18–33$ mag. Similar trends of this transition were found in the dense region of the Serpens South filament (Pillai et al. 2020) and in the Serpens Main dense cloud (Kwon et al. 2022), occurring at $A_V \sim 21$ mag and ~50 mag, respectively. This is likely because the B-fields are already dragged into the flow of very dense material during the gravitational contraction. Interestingly, we observe that the offset angle then progressively increases with $A_V \ge 40$ mag

<figure><img src="image_12.png" /></figure>

**Fig. 10.** Distribution of the offset angles between the orientations of B-fields and the ridge in ρ Oph-A, as a function of $A_V$ (circle symbols). The data in gray color represent the offset angle collected within the black circle “I” centered at the position of the Oph-S1 star with a radius of 60'' in Fig. 9. The data in blue and green additionally represent the offset angle collected within the two rings “II” and “III” in Fig. 9, respectively. The red diamonds represent the median value of the offset angles in each interval of $A_V=3$ mag. The number of red diamonds is reduced to improve the clarity of the figure. The red curve connecting all the red diamonds shows the trend of the offset angles vs. $A_V$ in ρ Oph-A. The cyan box indicates the range of $A_V=42–60$ mag where relative orientation between B-fields and the ridge changes from parallel to perpendicular. The gray horizontal dashed lines indicate perpendicular, 45°, and parallel relative orientations of B-fields with respect to the ridge.

in our region, indicating the relative orientation of B-field with respect to the ridge changes back from parallel to perpendicular (with the transition at $A_V \sim 42–60$ mag). We note, however, that this feature was not seen in the Serpens South filamentary cloud, whose maximum is $A_V \sim 60$ mag (see Fig. 3 in Pillai et al. 2020); however, it was seen in the Serpens Main cloud region, where $A_V \ge 170$ mag (see Fig. 6 in Kwon et al. 2022).

### 4.3. Magnetic field versus gas density

The B-field strength ($B_{tot}$ or $B$) can be expressed as a function of $n_H$ ($n_H = 2 n_{H_2}$): $B \propto n_H^k$. Theoretical models predict $k \sim 2/3$ in the weak B-field regions (Mestel 1966) and $k \lesssim 0.5$ in the strong B-field regions, where magnetic energy dominates over the gravity (Mouschovias & Ciolek 1999). Consequently, in the latter case, it suggests that $n_H$ would increase faster than $B$ if there is any collapse in the dense cloud.

Here, we show in Fig. 11 the relation $\log B - \log n_H$ toward the ρ Oph-A cloud, using B-fields measured with the DCF method in band D. We find that for densities $\le 10^{5.2} \text{ cm}^{-3}$, the B-fields distribution tends to be flattened. In a regime with densities higher than $10^{5.2} \text{ cm}^{-3}$ we start seeing the increasing trend of B-fields with respect to gas densities. Thus, we decided to test the $\log B - \log n_H$ relation at the densest part of the cloud only. We chose this region as the area within the third contour of Stokes-I at band D (see in the bottom right panel of Fig. 5), where Stokes-I $\ge 0.4$ Jy/arcsec$^2$. Within this region, it reveals a relatively weak correlation between the B-fields and densities (with the Pearson coefficient of the correlation is found

<footnote n="5">⁵ We convert $A_V$ from $N_{H_2}$ using the relation: $N_{H_2}/A_V = 0.94 \times 10^{21}$ molecules cm$^{-2}$ mag$^{-1}$ with a “standard” assumption of the total-to-selective extinction ratio, $R_V = A_V/E(B-V)$, of 3.1 (Frerking et al. 1982).</footnote>
<figure><img src="image_13.png" /><figcaption>Fig. 11. Correlation between log B and log n<sub>H</sub> in ρ Oph-A. Small circles represent B strength estimated with the DCF method at band D. Small circles in blue further represent data within the region where Stokes-I ≥ 0.4 Jy/arcsec<sup>2</sup>, corresponding to the densities log n<sub>H</sub> ≥ 5.2. The bigger circles in green represent median values of B binned in each interval of log n<sub>H</sub> of ~0.08 (see texts). Mean values of B estimated with ADF (red triangles), DMA (blue diamonds), and ST (orange squares) methods toward subregion (a) and (e) are overplotted. The solid green line represents the median value log B ≈ 2.81 at low densities. The solid blue line shows the power-law fit toward the binned data of B using the DCF method toward more denser region (covered by the orange shaded box), resulting k = 0.47 ± 0.18. The black dash-dotted line indicates the critical power-law B ∝ n<sub>H</sub><sup>0.5</sup> for strong B-fields, predicted from the theoretical models (Mouschovias & Ciolek 1999). The black dashed line indicates the relation found by Crutcher et al. (2010) using the Zeeman observational data.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Fig. 12. Image of the continuum emission at 154 µm toward Oph–A. B-field vectors obtained from SOFIA/HAWC+ at 154 µm are shown in magenta solid segments. The green dashed arrows indicate the radiation directions from the Oph-S1 star toward the cloud. The white star marks the position of the Oph-S1 star. The red curve illustrates the elongated shape of the ridge in the dense cloud. The gray contours are the H<sub>2</sub> column density, with contour levels at (1, 2, 5, 10, 20)×10<sup>22</sup> cm<sup>-2</sup>.</figcaption></figure>

r ~ 0.32, corresponding to a significance of ~5σ). We performed the power-law fit using the Python `curve_fit` routine (log B = k × log(n<sub>H</sub>) + b), using mean values of B estimated within every interval of log n<sub>H</sub>=0.08<sup>6</sup>. This gives a power-law index of k = 0.47 ± 0.18. In Fig. 11, we also overplot the results of B<sub>pos</sub> using other methods in subregions (a) and (e). We find that the fitted line well covers the data point representing the resulting B<sub>pos</sub> value using ADF method and shows good agreement with those estimated with ST and DMA in subregion (a) within and uncertainty of 1σ. The resulting k = 0.47 ± 0.18 is found consistent with theoretical predictions for strong B-fields (Mouschovias & Ciolek 1999). This is also consistent with β ≪ 1, representing strong B-fields in ρ Oph-A (see bottom panel in Fig. 8 and Sect. 3.6). We then performed the fit to B<sub>pos</sub><sup>ST</sup> values measured with the ST method from Sect. 3.3.2 and find k = 0.49 ± 0.27. We found a similar power-law index k for the dataset at band C, namely, 0.43 ± 0.12 and 0.39 ± 0.08, using result of B<sub>pos</sub><sup>DCF</sup> and B<sub>pos</sub><sup>ST</sup>, respectively (see Appendix C).

Our resulting k value is smaller than k ~ 0.65 found by Crutcher et al. (2010). However, their result was based on the mean values of B and n<sub>H</sub> in different molecular clouds, which contain significant observational uncertainty and has been questioned by later studies (for e.g., Tritsis et al. 2015; Li 2021).

## 4.4. Magnetic fields versus stellar feedback

High-mass stars irradiate, heat up, and ionize the surrounding gas and create the H II region (see e.g., Zinnecker & Yorke 2007). The ionized, hot gas also expands and interacts with the vicinity environment. The radiation pressure acting on both dust and gas (assuming dust is well-coupled to the gas in the cloud) can

drive the material away into the denser region. The lower density matter in the outer regions of ρ Oph-A is more exposed to the ionization front and pushed away faster than the higher density matter in the inner regions. This results in dust and gas being compressed and forming the dense shell structure at the edge of the cloud. Here, we investigate a possible impact of the radiation field from the high-mass star Oph-S1 on the ρ Oph-A cloud. The ρ Oph-A molecular cloud is directly exposed to the strong UV radiation field emitted by the Oph-S1 star from the east. Figure 12 shows the map of the continuum emission at 154 µm from SOFIA/HAWC+, superimposed by the white contours showing the H<sub>2</sub> column density.

We calculated the external radiation pressure of the Oph-S1 star acting on dust and gas in the densest part of the cloud, given by:

$$P_{\text{ext,rad}} = \frac{L_{\text{bol}}}{4\pi R^2 c}, \qquad (16)$$

where L<sub>bol</sub>=1100 L<sub>☉</sub> is the bolometric luminosity of the Oph-S1 star adopted from Lada & Wilking (1984) and R ≈ 0.06 pc is the projected distance from the cloud to the Oph-S1 star. We obtained P<sub>ext,rad</sub>/k<sub>B</sub> ≈ 2.4 × 10<sup>6</sup> K cm<sup>-3</sup>, which is in agreement with the result calculated from the ionizing photons emitted by the Oph-S1 star (~2×10<sup>6</sup> K cm<sup>-3</sup>, Pattle et al. 2015). Here, we omit the fact that the UV radiation from the Oph-S1 star can be absorbed by intervening dust in between the star and the ρ Oph-A region. Thus, the actual P<sub>ext,rad</sub> acting on dust and gas in ρ Oph-A is overestimated in the Eq. (16). To make a comparison with the external radiation pressure, we calculated the internal thermal pressure due to the hot gas in ρ Oph-A, P<sub>th</sub> = n<sub>H</sub><sub>2</sub>k<sub>B</sub>T<sub>gas</sub>. We considered the hot gas within the densest part of ρ Oph-A, where Stokes-I ≥ 0.4 Jy/arcsec<sup>2</sup> at 154 µm only, which corresponds to the area inside the third contour shown in

<footnote n="6"><sup>6</sup> The bin width is chosen by using the Freedman–Daconis rule.</footnote>
bottom right panel of Fig. 5. Within this region, we adopted $n_{H_2}$ in the range of $10^4–10^6$ cm$^{-3}$ and the average temperature of the hot gas $T_{gas} \approx T_{dust} \approx 20$ K, which yields $P_{th}/k_B$ in the range of $(0.2–20)\times10^6$ K cm$^{-3}$. With this assumption, we note that $P_{th}/k_B$ could be underestimated, since there is a gas heating process in the cloud dominated by the UV radiation from the Oph-S1 star and by the photoelectric emission from dust exposed to the UV radiation (i.e., $T_{gas} \ge T_{dust}$). The highest internal thermal pressure of $\sim 20\times10^6$ K cm$^{-3}$ is seen toward the densest part of the cloud, associated with the position of the SM1 core. This suggests that the internal thermal pressure is strong enough to help the material refrain from external pressure. We note that, however, in the vicinity region away from the center (e.g., the tails of the ridge shown in Fig. 12), the internal thermal pressure is weaker, at least one order of magnitude less than the external one. In this case, the radiation pressure from the Oph-S1 star on the east side is dominant and able to push away all the matter into the surrounding environment toward the west side. However, it is clear that the material in the low-density region still survives, despite owning a low internal thermal pressure (see Fig. 12). This indicates that there must be another mechanism at play, for example, the presence of strong B-fields of a few mG in the cloud. We roughly estimated the magnetic pressure ($P_{mag} = B^2/(8\pi)$) to be $10^7–10^8$ K cm$^{-3}$ for B-field strengths between 1–2 mG. This indicates that the magnetic pressure is at least two orders of magnitude higher than the external pressure of $\sim 10^6$ K cm$^{-3}$, sufficient to support the matter against the external pressure along with the gas thermal pressure.

We further investigated the interplay between B-fields and radiation directions on the evolution of the cloud. We determined the main radiation direction from the Oph-S1 star toward the densest region associated with the position of the starless core SM1 as the sky-projected direction from the Oph-S1 star and SM1 core. We found that the main radiation direction is $\sim 70^\circ$, measured east of north. The main B-field orientation toward the densest region of $\rho$ Oph-A of $\sim 74^\circ$ is almost parallel to the radiation field. The strong radiation pressure from the Oph-S1 star pushes all the material away and creates a dense shell toward the west side. To prevent dust and gas compression in the dense region, the B-fields exert a magnetic pressure perpendicular to the B-field orientation, seen perpendicular to the direction of the radiation pressure. This suggests that the magnetic pressure does not help the matter resist the external pressure, but tends to stretch the structure out along the northwest-southeast direction, resulting in an elongated shape of the dense shell. In the west side, the external thermal pressure is blocked by the gas thermal pressure due to the high-density material, resulting in the ridge shape in the cloud with very high densities. These results are consistent with results from 3D MHD simulation, which predicted that if the B-field orientation is parallel or along with the radiation direction, the cloud will tend to have a flattened shape Henney et al. (2009). However, in the northeast and southeast regions, the offset angles between the B-field vectors and radiation field are larger; namely, the B-field orientation becomes perpendicular to the radiation direction. In this case, the magnetic pressure and tension induced by B-fields run parallel with the radiation pressure, suggesting magnetic pressure can support the matter to confront the radiative compression and prevent the dust and gas from being swept away. This results in bending the material structure and forming a dense ridge with an umbrella-like shape. Interestingly, in the case study of $\rho$ Oph-A the ionization source is very close to the cloud ($<0.1$ pc) and the offset

position angle between the radiation and B-fields is not homogeneous; instead, it varies from parallel to orthogonal. This allows us to closely examine how the interplay between the B-fields and radiation fields could influence the material structure of the cloud.

## 5. Summary

In this paper, we present our analysis of the Ophiuchus dark cloud complex, using the SOFIA/HAWC+ polarimetric observations toward the $\rho$ Oph-A molecular cloud in the L1688 active star-forming region. Our main results are as follows:

- Using the dust polarization maps at 89 and 154 µm, we derived the map of B-field morphology in $\rho$ Oph-A. We found that B-fields are generally well-ordered and primarily perpendicular to the ridge of the cloud toward the higher density regime. We estimated $B_{pos}$ in the entire region of $\rho$ Oph-A using the DCF method, ranging from 0.2 to 2.3 mG at 89 and to 2.5 mG at 154 µm. We found that $B_{pos}$ is prominent toward the densest region of the cloud associated with the starless core SM1 and decreases toward the outskirts. This is the first time we have successfully estimated the distribution of the B-field strengths in $\rho$ Oph-A, which is essential to understanding how B-fields vary in the molecular cloud globally. Our derived B-field morphology and strengths well agree with previous measurements using the sub-mm dust polarization (Kwon et al. 2018). We also estimated the $B_{pos}$ strength using more advanced techniques, including the ADF, ST, and DMA, toward two subregions associated with the densest part of the cloud. While the $B_{pos}$ strengths estimated with the ADF method were in agreement with those from DCF, the ST and DMA methods provided lower $B_{pos}$ values of $\sim 4.5$ and 2–3 times, respectively;

- We studied the role of B-fields relative to gravity and turbulence by calculating the maps of the mass-to-flux ratio, Alfvénic Mach number, and plasma $\beta$ parameter toward $\rho$ Oph-A, using the inferred $B_{pos}^{DCF}$ from SOFIA/HAWC+ data at 154 µm. We found that toward the central dense parts the cloud is magnetically sub-critical and sub-Alfvénic where B-fields dominate over the gravity and turbulence. We obtained a result of $\beta \ll 1$ in $\rho$ Oph-A, suggesting the cloud is supported by the strong magnetic pressure and much greater than the thermal gas pressure. The lower $B_{pos}$ strengths obtained with the ST and DMA methods do not significantly change our conclusions on the roles of B-fields relative to gravity and turbulence on star formation. Our virial analysis in the entire cloud showed the virial parameter $\alpha_{vir} \sim 15$, suggesting that the cloud is gravitationally unbound and able to against the self-gravitational collapse. The cloud may be in the magneto hydrostatic equilibrium and may undergo gravitational collapse due to turbulent fluctuations on smaller scales. We should notice that our rough estimation of $\alpha_{vir}$ is a global parameter, applied for the entire cloud. In contrast, $\rho$ Oph-A is known as one of the densest parts of the active star-forming region L1688 hosting a few pre- and proto-stellar cores and embedded protostars. Higher spatial resolution observations toward the smaller scale region in the cloud will be necessary to investigate the dynamical state of the cloud in greater detail;

- We examined the role of B-fields in supporting the cloud against the radiation pressure from the nearby high-mass star Oph-S1. We found that at the center part of the cloud, the
magnetic pressure indeed does not play a role in helping the material to avoid the onset of radiation pressure. However, toward the northeast and southeast regions of the cloud, the magnetic pressure positively supports the cloud to prevent both the effects of radiative pressure and gas thermal pressure.

**Acknowledgements.** We thank our anonymous referee for their careful reading of our manuscript and valuable suggestions, which improved the quality of this work, particularly with regard to the realistic estimate of the gas volume density in the ρ Oph-A cloud. We thank Bhaswati Mookerjea for the valuable comments and discussion at the early stage of the project and for kindly sharing the final data cube of the HCO$^{+}$(4–3) transition line. We thank Fabio Santos for sharing the maps of dust temperatures and gas column density. We also thank Alex Lazarian for the constructive discussions concerning the DCF and DMA methods. N.L. thanks Nguyen Tat Thang for his valuable comments as part of the 3D structure of the ρ Oph-A. N.L. and A.K. acknowledge support from the First TEAM grant of the Foundation for Polish Science No. POIR.04.04.00-00-5D21/18-00 (PI: A. Karska). A.K. acknowledges also support from the Polish National Agency for Academic Exchange grant No. BPN/BEK/2021/1/00319/DEC/1. N.B.N. was funded by the Master, Ph.D. Scholarship Programme of Vingroup Innovation Foundation (VINIF), code VINIF.2023.TS.077. M.H. acknowledges support by the National Science Centre through the OPUS grant no. 2015/19/B/ST9/02959. N.L. also acknowledges support by the grant from the Simons Foundation to IFIRSE, ICISE (916424, PI: N.H.). P.N.D. and N.T.P. were funded by the Vietnam Academy of Science and Technology under project code THTETN.03/24-25. This work is based on observations made with the NASA/DLR Stratospheric Observatory for Infrared Astronomy (SOFIA). SOFIA is jointly operated by the Universities Space Research Association, Inc. (USRA), under NASA contract NNA17BF53C, and the Deutsches SOFIA Institut (DSI) under DLR contract 50 OK 2002 to the University of Stuttgart. This work made use of data from the JCMT/HARP and APEX. JCMT is operated by the East Asian Observatory on behalf of National Astronomical Observatory of Japan; the UK STFC under the auspices of grant number ST/N005856/1; Academia Sinica Institute of Astronomy and Astrophysics; the Korea Astronomy and Space Science Institute; the Operation, Maintenance and Upgrading Fund for Astronomical Telescopes and Facility Instruments, budgeted from the Ministry of Finance of China. This work used the facilities of the Canadian Astronomy Data Centre (CADC) operated by the National Research Council of Canada.

## References

Andre, P., Montmerle, T., Feigelson, E. D., Stine, P. C., & Klein, K. L. 1988, in *Activity in Cool Star Envelopes*, 143, eds. O. Havnes, B. R. Pettersen, J. H. M. M. Schmitt, & J. E. Solheim, 293

André, P., Belloche, A., Motte, F., & Peretto, N. 2007, A&A, 472, 519

Bertoldi, F., & McKee, C. F. 1992, ApJ, 395, 140

Borlaff, A. S., Lopez-Rodriguez, E., Beck, R., et al. 2023, ApJ, 952, 4

Brown, R. L., & Zuckerman, B. 1975, ApJ, 202, L125

Buckle, J. V., Hills, R. E., Smith, H., et al. 2009, MNRAS, 399, 1026

Chandrasekhar, S., & Fermi, E. 1953, ApJ, 118, 113

Chen, H. H.-H., Pineda, J. E., Goodman, A. A., et al. 2019, ApJ, 877, 93

Chuss, D. T., Andersson, B. G., Bally, J., et al. 2019, ApJ, 872, 187

Connelley, M. S., & Greene, T. P. 2010, AJ, 140, 1214

Coudé, S., Bastien, P., Houde, M., et al. 2019, ApJ, 877, 88

Crutcher, R. M. 2004, in *The Magnetized Interstellar Medium*, eds. B. Uyaniker, W. Reich, & R. Wielebinski, 123

Crutcher, R. M. 2012, ARA&A, 50, 29

Crutcher, R. M., Wandelt, B., Heiles, C., Falgarone, E., & Troland, T. H. 2010, ApJ, 725, 466

Davis, L. 1951, *Phys. Rev.*, 81, 890

de Zeeuw, P. T., Hoogerwerf, R., de Bruijne, J. H. J., Brown, A. G. A., & Blaauw, A. 1999, AJ, 117, 354

Draine, B. T., & Weingartner, J. C. 1997, ApJ, 480, 633

Enoch, M. L., Evans, Neal J., I., Sargent, A. I., & Glenn, J. 2009, ApJ, 692, 973

Esplin, T. L., & Luhman, K. L. 2020, AJ, 159, 282

Evans, Neal J., I., Dunham, M. M., Jørgensen, J. K., et al. 2009, ApJS, 181, 321

Frerking, M. A., Langer, W. D., & Wilson, R. W. 1982, ApJ, 262, 590

Gordon, M. S., Lopez-Rodriguez, E., Andersson, B. G., et al. 2018, arXiv e-prints [arXiv:1811.03100]

Greve, T. R., Papadopoulos, P. P., Gao, Y., & Radford, S. J. E. 2009, ApJ, 692, 1432

Guerra, J. A., Chuss, D. T., Dowell, C. D., et al. 2021, ApJ, 908, 98

Güsten, R., Nyman, L. Å., Schilke, P., et al. 2006, A&A, 454, L13

Hall, J. S. 1949, *Science*, 109, 166

Hamaguchi, K., Corcoran, M. F., & Imanishi, K. 2003, PASJ, 55, 981

Harper, D. A., Runyan, M. C., Dowell, C. D., et al. 2018, J. Astron. Instrum., 7, 1840008

Hennebelle, P., & Falgarone, E. 2012, A&A Rev., 20, 55

Hennebelle, P., & Inutsuka, S.-i. 2019, Front. Astron. Space Sci., 6

Henney, W. J., Arthur, S. J., de Colle, F., & Mellema, G. 2009, MNRAS, 398, 157

Hildebrand, R. H., Kirby, L., Dotson, J. L., Houde, M., & Vaillancourt, J. E. 2009, ApJ, 696, 567

Hiltner, W. A. 1949, Nature, 163, 283

Hoang, T., & Lazarian, A. 2016, ApJ, 831, 159

Hoang, T., Tram, L. N., Minh Phan, V. H., et al. 2022, AJ, 164, 248

Holland, W. S., Bintley, D., Chapin, E. L., et al. 2013, MNRAS, 430, 2513

Hwang, J., Kim, J., Pattle, K., et al. 2021, ApJ, 913, 85

Hwang, J., Pattle, K., Parsons, H., Go, M., & Kim, J. 2023, AJ, 165, 198

Kandori, R., Kusakabe, N., Tamura, M., et al. 2006, in *Ground-based and Airborne Instrumentation for Astronomy*, 6269, eds. I. S. McLean & M. Iye, International Society for Optics and Photonics (SPIE), 626951

Kauffmann, J., Bertoldi, F., Bourke, T. L., Evans, N. J., I., & Lee, C. W. 2008, A&A, 487, 993

Kennicutt, R. C., & Evans, N. J. 2012, ARA&A, 50, 531

Koch, E. W., & Rosolowsky, E. W. 2015, MNRAS, 452, 3435

Krumholz, M. R., Bate, M. R., Arce, H. G., et al. 2014, in Protostars and Planets VI, eds. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 243

Kwon, J., Tamura, M., Hough, J. H., et al. 2015, ApJS, 220, 17

Kwon, J., Doi, Y., Tamura, M., et al. 2018, ApJ, 859, 4

Kwon, W., Pattle, K., Sadavoy, S., et al. 2022, ApJ, 926, 163

Lada, C. J., & Wilking, B. A. 1984, ApJ, 287, 610

Lada, C. J., & Lada, E. A. 2003, ARA&A, 41, 57

Ladjelate, B., André, P., Könyves, V., et al. 2020, A&A, 638, A74

Lazarian, A., & Vishniac, E. T. 1999, ApJ, 517, 700

Lazarian, A., & Hoang, T. 2007, MNRAS, 378, 910

Lazarian, A., Yuen, K. H., & Pogosyan, D. 2022, ApJ, 935, 77

Li, H.-B. 2021, *Galaxies*, 9, 41

Liseau, R., Justtanont, K., & Tielens, A. G. G. M. 2006, A&A, 446, 561

Liseau, R., Larsson, B., Lunttila, T., et al. 2015, A&A, 578, A131

Liu, J., Zhang, Q., Commerçon, B., et al. 2021, ApJ, 919, 79

Liu, J., Zhang, Q., & Qiu, K. 2022, Front. Astron. Space Sci., 9, 943556

Lopez-Rodriguez, E., Beck, R., Clark, S. E., et al. 2021, ApJ, 923, 150

Lopez-Rodriguez, E., Mao, S. A., Beck, R., et al. 2022, ApJ, 936, 92

Loren, R. B., Wootten, A., & Wilking, B. A. 1990, ApJ, 365, 269

Mackey, J., & Lim, A. J. 2011a, MNRAS, 412, 2079

Mackey, J., & Lim, A. J. 2011b, MNRAS, 412, 2079

Matthews, B. C., McPhee, C. A., Fissel, L. M., & Curran, R. L. 2009, ApJS, 182, 143

McKee, C. F., & Ostriker, E. C. 2007, ARA&A, 45, 565

Mestel, L. 1966, MNRAS, 133, 265

Mestel, L., & Spitzer, L. J. 1956, MNRAS, 116, 503

Mookerjea, B., Sandell, G., Vacca, W., Chambers, E., & Güsten, R. 2018, A&A, 616, A31

Motte, F., Andre, P., & Neri, R. 1998, A&A, 336, 150

Mouschovias, T. 2001, in *Astronomical Society of the Pacific Conference Series*, 248, Magnetic Fields Across the Hertzsprung-Russell Diagram, eds. G. Mathys, S. K. Solanki, & D. T. Wickramasinghe, 515

Mouschovias, T. C., & Ciolek, G. E. 1999, in *NATO Advanced Study Institute (ASI) Series C*, 540, The Origin of Stars and Planetary Systems, eds. C. J. Lada, & N. D. Kylafis, 305

Mouschovias, T. C., Tassis, K., & Kunz, M. W. 2006, ApJ, 646, 1043

Nakano, T., & Nakamura, T. 1978, PASJ, 30, 671

Ngoc, N. B., Diep, P. N., Parsons, H., et al. 2021, ApJ, 908, 10

Ngoc, N. B., Diep, P. N., Hoang, T., et al. 2023, ApJ, 953, 66

Ortiz-León, G. N., Loinard, L., Kounkel, M. A., et al. 2017, ApJ, 834, 141

Ostriker, E. C., Stone, J. M., & Gammie, C. F. 2001, ApJ, 546, 980

Pabst, C., Higgins, R., Goicoechea, J. R., et al. 2019, Nature, 565, 618

Pabst, C. H. M., Goicoechea, J. R., Teissier, D., et al. 2020, A&A, 639, A2

Pattle, K., & Fissel, L. 2019, Front. Astron. Space Sci., 6, 15

Pattle, K., Ward-Thompson, D., Kirk, J. M., et al. 2015, MNRAS, 450, 1094

Pattle, K., Ward-Thompson, D., Berry, D., et al. 2017, ApJ, 846, 122

Pattle, K., Ward-Thompson, D., Hasegawa, T., et al. 2018, ApJ, 860, L6

Pattle, K., Lai, S.-P., Di Francesco, J., et al. 2021, ApJ, 907, 88

Pattle, K., Lai, S.-P., Sadavoy, S., et al. 2022, MNRAS, 515, 1026

Pattle, K., Fissel, L., Tahani, M., Liu, T., & Ntormousi, E. 2023, in *Astronomical Society of the Pacific Conference Series*, 534, Protostars and Planets VII, eds. S. Inutsuka, Y. Aikawa, T. Muto, K. Tomida, & M. Tamura, 193

Pillai, T., Kauffmann, J., Wyrowski, F., et al. 2011, A&A, 530, A118

Pillai, T. G. S., Clemens, D. P., Reissl, S., et al. 2020, Nat. Astron., 4, 1195

Planck Collaboration Int. XXXV. (Ade, P. A. R., et al.) 2016, A&A, 586, A138

Plummer, H. C. 1911, MNRAS, 71, 460

Poglitsch, A., Waelkens, C., Geis, N., et al. 2010, A&A, 518, L2
Poidevin, F., Falceta-Gonçalves, D., Kowal, G., de Gouveia Dal Pino, E., & Mário Magalhães, A. 2013, ApJ, 777, 112

Purcell, E. M. 1979, ApJ, 231, 404

Santos, F. P., Chuss, D. T., Dowell, C. D., et al. 2019, ApJ, 882, 113

Santos-Lima, R., Dal Pino, E. d. G., & Lazarian, A. 2012, ApJ, 747, 21

Shu, F. H., Adams, F. C., & Lizano, S. 1987, ARA&A, 25, 23

Skalidis, R., & Tassis, K. 2021, A&A, 647, A186

Skalidis, R., Sternberg, J., Beattie, J. R., Pavlidou, V., & Tassis, K. 2021, A&A, 656, A118

Soam, A., Maheswar, G., Lee, C. W., Neha, S., & Kim, K.-T. 2018, MNRAS, 476, 4782

Soler, J. D. 2019, A&A, 629, A96

Soler, J. D., Ade, P. A. R., Angilè, F. E., et al. 2017, A&A, 603, A64

Tram, L. N., Hoang, T., Lee, H., et al. 2021, ApJ, 906, 115

Tram, L. N., Bonne, L., Hu, Y., et al. 2023, ApJ, 946, 8

Tritsis, A., Panopoulou, G. V., Mouschovias, T. C., Tassis, K., & Pavlidou, V. 2015, MNRAS, 451, 4384

Vaillancourt, J. E. 2006, PASP, 118, 1340

Ward-Thompson, D., Pattle, K., Bastien, P., et al. 2017, ApJ, 842, 66

Ward-Thompson, D., Karoly, J., Pattle, K., et al. 2023, ApJ, 946, 62

White, G. J., Drabek-Maunder, E., Rosolowsky, E., et al. 2015, MNRAS, 447, 1996

Wilking, B. A., Gagné, M., & Allen, L. E. 2008, in Handbook of Star Forming Regions, Volume II, 5, ed. B. Reipurth, 351

Zinnecker, H., & Yorke, H. W. 2007, ARA&A, 45, 481

Zucker, C., & Chen, H. H.-H. 2018, ApJ, 864, 152
# Appendix A: Histogram of SOFIA/HAWC+ dataset

Figures A.1 and A.2 show histograms of Stokes-I, $\sigma_I$, SNR$_I$, p, $\sigma_p$, and SNR$_p$ of the SOFIA/HAWC+ data in band C and band D. Table A.1 shows the summary of the mean and RMS of the distribution of Stokes-I, p, their uncertainties and SNRs.

<figure><img src="image_15.png" /></figure>

Fig. A.1: Histograms of Stokes-I (top panels), $\sigma_I$ (middle panels), and SNR$_I$ (bottom panels) of the SOFIA/HAWC+ data. Left panels represent data in band C at 89 µm and right panels represent data in band D at 154 µm. The blue vertical dashed lines in each panel represent the mean value of the distribution. The red vertical dashed lines show the position of mean value plus 1$\sigma$, 2$\sigma$, and 3$\sigma$ of the distribution from left to right, respectively.

<figure><img src="image_16.png" /></figure>

Fig. A.2: Histograms of p (top panels), $\sigma_p$ (middle panels), and SNR$_p$ (bottom panels) of the SOFIA/HAWC+ data. Left panels represent data in band C at 89 µm and right panels represent data in band D at 154 µm. The blue vertical dashed lines in each panel represent the mean value of the distribution. The red vertical dashed lines show the position of mean value plus 1$\sigma$, 2$\sigma$, and 3$\sigma$ of the distribution from left to right, respectively.

# Appendix B: Magnetic field orientations and strengths at different wavelengths

In Fig. B.1, we compare the B-fields vectors obtained from HAWC+ data at 89 and 154 µm (this study), and from POL-2 data at 850 µm (Kwon et al. 2018). Vectors show the mean B-field position angles at the subregions in $\rho$ Oph-A.

Figure B.2 shows the histogram distribution of the estimated $B_{\text{pos}}^{\text{DCF}}$ obtained with HAWC+ data at 89 and 154 µm. Figure B.3 shows the maps of B-field strength estimated using the ST method ($B_{\text{pos}}^{\text{ST}}$) at 89 and 154 µm.

<figure><img src="image_17.png" /></figure>

Fig. B.1: Same as the right panel in Fig. 2 but with B-field vectors from POL-2 data added (red dash-dotted vectors, Kwon et al. 2018). These vectors represent the mean B-field position angles in subregions (see texts). White and black dashed vectors show the corresponding mean B-field position angles at 89 and 154 µm.

# Appendix C: Velocity dispersion traced by other molecular lines

To test whether the FWHM map of HCO$^+$ (4–3) line is a good choice, we further use the N$_2$H$^+$ (3–2) molecular data observed
<table><caption>Table A.1: Mean and RMS of Stokes-I, *p*, and their SNR in bands C and D.</caption><thead><tr><th rowspan="2">Wavelength</th><th colspan="2">I (Jy/arcsec<sup>2</sup>)</th><th colspan="2">σ<sub>I</sub> (Jy/arcsec<sup>2</sup>)</th><th colspan="2">SNR<sub>I</sub></th><th colspan="2">p (%)</th><th colspan="2">σ<sub>P</sub> (%)</th><th colspan="2">SNR<sub>P</sub></th></tr><tr><th>Mean</th><th>RMS</th><th>Mean</th><th>RMS</th><th>Mean</th><th>RMS</th><th>Mean</th><th>RMS</th><th>Mean</th><th>RMS</th><th>Mean</th><th>RMS</th></tr></thead><tbody><tr><td>89 µm</td><td>0.163</td><td>0.166</td><td>0.005</td><td>0.004</td><td>44.8</td><td>52.9</td><td>8.3</td><td>9.2</td><td>6.3</td><td>7.8</td><td>3.2</td><td>3.5</td></tr><tr><td>154 µm</td><td>0.143</td><td>0.142</td><td>0.001</td><td>0.001</td><td>177.5</td><td>214.1</td><td>5.7</td><td>4.1</td><td>1.2</td><td>1.5</td><td>10.3</td><td>9.9</td></tr></tbody></table>

<figure><img src="image_18.png" /><figcaption>Fig. B.2: Histogram distribution of the estimated $B_{\text{pos}}$ using the DCF method, obtained with HAWC+ data at 89 µm (in blue color with hatches) and 154 (in orange color) µm.</figcaption></figure>

by APEX to estimate the total velocity dispersion, $\sigma_v$. The majority of the $\rho$ Oph-A cloud exhibits the $N_2H^+(3-2)$ emission showing 3 single Gaussian profiles; the middle peak is much stronger compared with the other two components. This strongest hyperfine line indicates the central velocity which is in the range of 2.54–3.97 km s$^{-1}$, with an average velocity of 3.42 km s$^{-1}$, in agreement with the systemic velocity of $\rho$ Oph-A of ~3.44 km s$^{-1}$ estimated by André et al. (2007). To each pixel we fit the spectrum with the sum of 3 single Gaussian profiles and then use the line width of the primary peak as the $\sigma_v$ toward the computing pixel. We construct the map of the FWHM velocity dispersion, following the procedure presented in Sect. 3.2.2 for the case of the HCO$^+(4-3)$ molecular line.

Figure C.1 shows the map of $\Delta v_{\text{NT}}$ using the $N_2H^+(3-2)$ molecular line toward $\rho$ Oph-A. The $\Delta v_{\text{NT}}$ values ranges from 0.18 to 1.73 km s$^{-1}$ with a median value of 0.95 km s$^{-1}$, in good agreement with those obtained by HCO$^+(4-3)$ line. We note, however, that the spatial extent of the HCO$^+(4-3)$ emission line is larger than that of the $N_2H^+(3-2)$. Thus, we use the map of FWHM of the HCO$^+(4-3)$ line to calculate the $B_{\text{pos}}^{\text{DCF}}$ strengths in $\rho$ Oph-A. Similarly, White et al. (2015) reported a mean FWHM velocity dispersion of 1.5 km s$^{-1}$ using the C$^{18}$O(3–2) data from JCMT/HARP in overall entire $\rho$ Ophiuchus cloud. Using the $N_2H^+(1-0)$ data from IRAM, André et al. (2007) estimated an average value of the 1D nonthermal velocity dispersion toward $\rho$ Oph-A, $\sigma_{\text{NT}}=0.2$ km s$^{-1}$, corresponding to $\Delta v_{\text{NT}} \sim 0.47$ km s$^{-1}$. Observations of the $\text{NH}_3$ hyperfine emission from the Green Bank Telescope also returned a similar 1D velocity dispersion with the average value of 0.3 km s$^{-1}$, corresponding to an average of FWHM of 0.65 km s$^{-1}$ (Chen et al. 2019).
<figure><img src="image_19.png" /></figure>

Fig. B.3: Maps of the B-field strength projected in the plane-of-the-sky in $\rho$ Oph-A obtained by using the ST method, $B_{\text{pos}}^{\text{ST}}$, in band C (left panel) and band D (right panel). Other notations are similar as in Fig. 6.

<figure><img src="image_20.png" /></figure>

Fig. C.1: Distribution of the non-thermal velocity FWHM of the $N_2H^+$ (3–2) molecular line toward $\rho$ Oph-A. White star marks the position of the Oph-S1 star. White contours show the continuum emission (Stokes-I) at 154 µm, with contours levels at 0.12, 0.2, 0.4, 0.6, and 0.8 Jy/arcsec$^2$.

## Appendix D: Mass-to-flux ratio, Alfvénic Mach number, beta plasma, and $B-n_H$ relation results obtained by using $B_{\text{pos}}^{\text{ST}}$

Figure D.1 shows the maps of mass-to-flux ratio (top panel), Alfvénic Mach number (middle panel), and beta plasma (bottom panel) toward $\rho$ Oph-A using $B_{\text{pos}}^{\text{DCF}}$ in band C. Figure D.2 shows maps of mass-to-flux ratio (top), Alfvénic Mach number (middle), and beta plasma (bottom) toward $\rho$ Oph-A using $B_{\text{pos}}^{\text{ST}}$ in band C (left panels) and D (right panels).

Figure D.3 shows the relation between $\log B$ and $\log n_H$ using $B_{\text{pos}}^{\text{ST}}$ in band D.
<figure><img src="image_21.png" /><figcaption>Fig. D.1: Similar to Fig. 8 but using $B_{\text{pos}}^{\text{DCF}}$ in band C.</figcaption></figure>
<figure><img src="image_22.png" /><figcaption>Fig. D.2: Similar to Fig. 8 but using $B_{pos}^{ST}$ in band C (left panels) and D (right panels).</figcaption></figure>
<figure><img src="image_23.png" /></figure>

Fig. D.3: Similar to Fig. 11 but using B-field strength estimated with the ST method in band D. Small squares represent *B* strength estimated with the ST method at band D. Blue small squares further show the data at the denser region (log *n*<sub>H</sub> ≥ 5.2) with Stokes-*I* ≥ 0.4 Jy/arcsec². Green bigger squares represent median value of log *B* binned in each interval of log *n*<sub>H</sub> of ~0.08 in the denser region. Mean values of log *B* estimated with ADF (red triangles) and DMA (blue diamonds) methods toward subregions (a) and (e) are overplotted. The solid green line shows the median data at lower densities region where log *n*<sub>H</sub> ≥5.2. The solid blue line shows the power-law fit toward the binned data of *B* using the DCF method toward denser region (covered by the orange shaded box, resulting *k* = 0.49 ± 0.27. The black dash-dotted line indicates the critical power-law *B* ∝ *n*<sub>H</sub><sup>0.5</sup> for strong B-fields, predicted from the theoretical models (Mouschovias & Ciolek 1999). The black dashed line indicates the relation found by Crutcher et al. (2010) using the Zeeman observational data.

<figure><img src="image_24.png" /></figure>

Fig. E.1: Schematic view of the ridge in ρ Oph-A. The background color map shows the Stokes-*I* intensity at 154 µm. The thick red line indicates the main shape of the ridge. Thin red lines show the lines perpendicular to the tangent line at each position on the ridge.

**Appendix E: Main structure of the ridge
identification**

We identified the elongated shape of the ridge in the densest
region of ρ Oph-A using the Python package RadFil. We first
applied a threshold of Stokes-I ≥0.5 Jy/arcsec² in the map of
Stokes-I at 154 µm to build a mask for the main shape of the
densest region (the ridge). The spine of the ridge was identi-
fied by using the Python package FilFinder<fnref n="7" /> from Koch &
Rosolowsky (2015). We then smoothed the ridge spine and define
a tangent line at every position on the ridge, where the ori-
entation of the ridge (θ<sub>ridge</sub>) at every position was identified.
Figure E.1 shows the schematic view of the ridge and the per-
pendicular cuts at each position on the ridge. The position angle
of the tangent line at each position located on the ridge is defined
as the orientation of the ridge at that position, measured east of
north, and ranges from 0° to 180°.

<footnote n="7">7 https://github.com/e-koch/FilFinder</footnote>