Main

In non-masing interstellar medium, only H I, OH and CN have successfully produced systematic Zeeman measurements. The comprehensive Zeeman surveys15 indicate that the magnetic fields in the diffuse cold neutral medium (CNM) probed by H I do not scale significantly with density, whereas above a critical break-point density of approximatey 300 cm−3, the magnetic fields in dense cores probed by OH tend to increase with density. However, owing to the gap in densities between the H I (about 40 cm−3) and OH (\(\gtrsim \)103 cm−3) Zeeman measurements, the field transition around the critical density of about 300 cm−3 (where the dependence of the field strength on the density changes behaviour) remains a controversial topic and could have crucial implications for star formation16,17,18,19. Recently, a CCS Zeeman detection20 shed light on regions denser than that probed by OH. A Zeeman probe that is sensitive for a wide range of densities, particularly the low-density molecular envelope, is highly desirable and could help to distinguish different core-formation scenarios.

We developed the so-called H I narrow self-absorption (HINSA) technique to provide a probe of the transition from H I to H2 (refs. 4,5). HINSA traces cold atomic hydrogen well mixed with H2, which provides the necessary cooling, not available in the CNM, of H I through collision. Close to the steady state between H2 formation and destruction, the HINSA strength is independent of the gas density5 and thus capable of probing the transition around the critical density. Although the Zeeman effect of the H I self-absorption feature has been reported21,22, the broad line widths of the absorption components are mostly associated with diffuse atomic gas rather than dense molecular gas. Considering that HINSA typically has much higher brightness temperatures than most molecular lines, is impervious to depletion23 and can be detected in a wide range of H2 densities, HINSA is a promising Zeeman probe for molecular gas.

The HINSA feature in L1544 has a strong absorption dip and a nearly thermalized narrow line width at a temperature lower than 15 K (ref. 4). The non-thermal line width and centroid velocity of the HINSA are very close to those of the emission lines of OH, 13CO and C18O molecules, and their column densities are well correlated, suggesting that a significant fraction of the atomic hydrogen is located in the cold, well-shielded portions of L15445. We thus assume that the column density sampled by the HINSA can be approximated by that obtained from dust, despite the substantially larger apparent area covered by HINSA (Fig. 1a). The previous OH Zeeman detection with Arecibo24 towards the L1544 centre resulted in a field strength of Blos = +10.8 ± 1.7 μG, where Blos is the magnetic field component along the line of sight, with a positive sign representing the field pointing away from the observer. In contrast, the OH Zeeman observations of the Green Bank Telescope (GBT) towards four envelope locations 6.0′ (0.24 pc) from the centre yielded a marginal detection of Blos = +2 ± 3 μG (ref. 16), leaving the structure of the envelope field undetermined.

Fig. 1: L1544 core and illustration of the structure of interstellar medium from the CNM to the core.
figure 1

a, A composite of Digitized Sky Survey 2 (DSS2) images of L1544 with the i band in red, the r band in green and the b band in blue overlaid with the HINSA and H2 column density maps. The white dashed contours are 30%, 50%, 70% and 90% of the peak HINSA column density, and the orange contours are 2 × 1021 cm−2, 4 × 1021 cm−2, 6 × 1021 cm−2, 8 × 1021 cm−2 and 10 × 1021 cm−2 for the H2 column density. The red, green and cyan circles mark the locations and beam sizes of the FAST, Arecibo and GBT Zeeman observations, respectively. b, A composite of three 0.5 km s−1 velocity slices of the Arecibo GALFA-H I images at 6.2 km s−1, 6.7 km s−1 and 7.3 km s−1 local standard of rest (LSR) velocities in blue, green and red. The dashed rectangle shows the region of a. The two absorption dots represent the locations of quasars 3C132 and 3C133. c, Schematic view of the CNMs, the molecular envelope and the L1544 core.

With the Five-hundred-meter Aperture Spherical radio Telescope (FAST)25, we detected Zeeman splittings in a 2.9′ beam (0.12 pc) towards the HINSA column density peak, 3.6′ (0.15 pc) away from the L1544 centre (Fig. 1). The spectra of the Stokes I(v) and V(v) parameters (where v denotes velocity) are shown in Fig. 2. The I(v) spectrum contains the H I emission of the CNM and the warm neutral medium (WNM) clouds in the direction towards the Taurus complex and a HINSA feature at the centroid velocity of L1544. Figure 2a shows our decomposition of I(v) into a foreground HINSA component, a background WNM component, and three CNM components between the HINSA and WNM. Our fitted parameters of the HINSA component are in good agreement with the previous HINSA observations4,5, and our parameters of the CNM and WNM components are similar to the Arecibo results towards quasars around L154426.

Fig. 2: The Stokes I(v) and V(v) spectra at 21-cm wavelength towards the HINSA column density peak.
figure 2

a, The black profile represents the I(v) spectrum. The red profile represents the absorption from the foreground HINSA component. The blue dashed and dotted profiles represent the emission of the CNM and WNM components, respectively. The CNM and WNM profiles include the absorption from the CNM components that lie in front but do not include the absorption from the HINSA. The black dashed profile represents the sum of the absorption and emission profiles. b, The black profile represents the V(v) spectrum. The black dashed profile represents the sum of the Zeeman splitting profiles of the five components. The red profile represents the Zeeman splitting profile with Blos = +3.8 μG of the HINSA component.

The V(v) spectrum shows features of classic ‘S curve’ patterns proportional to the first derivatives of I(v) for the HINSA, CNM and WNM components, as expected for Zeeman splittings. The Zeeman splitting profile of HINSA has a maximum at high velocity and a minimum at low velocity, opposite to the Zeeman splitting profile of CNM1, the closest CNM component at a velocity similar to L1544, which shows positive V at low velocity and negative V at high velocity. From our least-squares fits to V(v), Fig. 2b shows the Zeeman splitting of the HINSA and the total Zeeman profile of the five components, and Fig. 3 shows the individual Zeeman splittings and Blos of the components. The HINSA Zeeman effect gives Blos = +3.8 ± 0.3 μG, and the H I Zeeman effect of CNM1 gives Blos = +4.0 ± 1.1 μG. The magnetic field strengths of HINSA and CNM1 are consistent with the results of Blos = +5.8 ± 1.1 μG and Blos = +4.2 ± 1.0 μG obtained from the Zeeman observations towards quasars 3C133 and 3C132, probing the magnetic fields of CNM1 at distances of 17.7′ (0.72 pc) and 174.5′ (7.1 pc) from L1544, respectively27. For the second and third CNM components (CNM2 and CNM3) along the line of sight, our results of Blos,CNM2 = −7.6 ± 1.0 μG and Blos,CNM3 = +2.9 ± 0.4 μG are also consistent with the results of Blos,CNM2 = −9.6 ± 6.3 μG and Blos,CNM3 = −0.3 ± 1.7 μG towards quasar 3C13327.

Fig. 3: The individual V(v) profiles for the HINSA, CNM and WNM components.
figure 3

In each panel, the red profile represents the fitted Zeeman profile of the component and the black profile represents the observed V(v) subtracted by the fitted Zeeman profiles of the other four components. The CNM and WNM Zeeman profiles include the absorption from the CNM components that lie in front but do not include the absorption from the HINSA. The sum of the red profiles of the five components is the black dashed profile in Fig. 2b.

Comparing the Zeeman observations of HINSA, OH and H I tracing the CNM1 and the molecular envelope of L1544, it is clear that the magnetic fields at distances of 0.15 pc, 0.24 pc, 0.72 pc and 7.1 pc from the centre all have the same direction of Blos and consistent strengths roughly within the 1σ. This finding is in agreement with the conclusion of a median value of 6 μG in absolute total strength in H I clouds inferred from comprehensive Zeeman surveys15. The HINSA Zeeman effect thus provides a connection between the magnetic fields from H I clouds to molecular clouds. The H I emission components (CNM1, CNM2 and CNM3) and the H I absorption towards 3C133 and 3C132 trace the CNM with a kinematic temperature of about 100 K (ref. 27) and a number density of about 40 cm−3 (ref. 15), whereas the HINSA and OH observations trace the envelope of about 10–15 K (ref. 4) and about 103 cm−3 (refs. 5,15). Despite the one to two orders of magnitude change in both temperature and density in the phase transition from the atomic CNM to the molecular envelope, the Zeeman observations reveal a magnetic field that is coherent in both direction and strength across multi-scales and multi-phases of the interstellar medium. To constrain the uniformity of the coherent magnetic field, our likelihood analysis of the HINSA, OH and H I Zeeman measurements suggests a Gaussian distribution of the Blos with a mean strength of B0 = +4.1 ± 1.6 μG and an intrinsic spread of \({\sigma }_{0}={1.2}_{-0.6}^{+1.2}\,{\rm{\mu }}{\rm{G}}\), a significantly better constraint than the previous estimation of \({B}_{0}=+{4}_{-8}^{+10}\,{\rm{\mu }}{\rm{G}}\) based on only the OH results18.

It is well known that the progenitor of molecular gas, the atomic CNM, is strongly magnetized, as measured by the dimensionless mass-to-magnetic-flux ratio λ in units of the critical value of 2πG1/2 (λ = 7.6 × 10−21[NH2 (cm−2)][Btot (μG)]−1, where NH2 is the column density of H2 gas and Btot is the total magnetic field strength1), which is well below unity (that is, magnetically subcritical)28. However, the immediate progenitors of stars, the prestellar cores of molecular clouds such as the L1554 core, are observed to be magnetically supercritical (λ > 1)24, which is required for the self-gravity to overwhelm the magnetic support and form stars through gravitational collapse. When and how the transition from the magnetically subcritical CNM that is incapable of forming stars through direct gravitational collapse to the supercritical star-forming cores occurs is a central unresolved question in star formation.

Our HINSA Zeeman observations can be used to address this question. Using the physical parameters of the clouds (Table 1) and the statistically most probable value of Btot = 2Blos, the λ of CNM1 is about 0.10–0.18, consistent with previous results28. The λ of the envelope and core of L1544 core is 2.5–3.5, which is well above unity, indicating that the transition to magnetic supercriticality has already occurred. We further consider the relative values of λ between CNM1 and L1544 to avoid the geometrical correction from Blos to Btot (ref. 16), assuming that the inclination angles of the magnetic fields in the L1544 core and envelope are similar. Therefore, the molecular envelope of the L1544 core traced by HINSA is at least 13 times less magnetized relative to its mass compared with its ambient CNM. This is different from the ‘classic’ theory of low-mass star formation, which envisions the transition from magnetic subcriticality to supercriticality occurring as the supercritical core forms out of the magnetically supported (subcritical) envelope13,14. Our results suggest that the transition from magnetic subcriticality to supercriticality occurs earlier, during the formation of the molecular envelope, favouring the more rapidly evolving scenario of core formation and evolution for L15448 over the slower, magnetically retarded scenario9. In other words, by the time that the molecular envelope is formed, the problem of excessive magnetic flux as a fundamental obstacle to gravitational collapse and star formation is already resolved. This early reduction of flux relative to mass is unlikely owing to the ‘classical’ scenario where gravity drives neutrals through ions (and the magnetic field tied to them) in a process called ‘ambipolar diffusion’ because the CNM is not self-gravitating. The coherent magnetic fields reviewed here provide a new specific question on how to create supercritical dense cores such as L1544 from subcritical clouds. Plausible scenarios include mass accumulation along field lines29 and (turbulence enhanced) magnetic reconnection30, although whether such scenarios can reproduce the distributions of gas and magnetic field observed in the L1544 region remains to be seen. In any case, the already magnetically supercritical envelope can in principle go on to form dense cores and stars without having to further reduce its magnetic flux relative to the mass.

Table 1 Physical parameters of the clouds

Methods

Data reduction

The FAST Zeeman observations towards the HINSA column density peak in L1544 were carried out on five days between August and November 2019 with a total integration time of 7.6 h. The HINSA spectra were obtained with the central beam of the L-band 19-beam receiver31. The central beam has an average system temperature of 24 K, a main beam efficiency of 0.63 and a main beam diameter at the half-power point of 2.9′ with a pointing accuracy of 7.9″. The 19-beam receiver had orthogonal linear polarization feeds followed by a temperature-stabilized noise injection system and low noise amplifiers to produce the X and Y signals of the two polarization paths. The XX, YY, XY and YX correlations of the signals then were simultaneously recorded using the ROACH backend with 65,536 spectral channels in each polarization. The spectral bandwidth was 32.75 MHz centred at the frequency of the H I 21-cm line for a channel spacing of 500 Hz, and the V(v) spectrum presented in this work was Hanning smoothed, which produced a spectral resolution of 0.21 km s−1.

The data reduction, including gain and phase calibrations of the two polarization paths, bandpass calibrations of the four correlated spectra and polarization calibrations to generate the Stokes I, Q, U and V spectra, was carried out using the IDL RHSTK package written by C. Heiles and T. Robishaw, which is widely used for Arecibo and GBT polarization data. The 19-beam receiver is rotatable from −80° to +80° with respect to the line of equatorial latitude. The polarization calibrations used drifting scans of the continuum source 3C286 at rotation angles of −60°, −30°, 0°, 30° and 60° over 1.5 h surrounding its transit. The details of the polarization calibration procedure are provided in ref. 32. We performed polarization calibrations once per month during the observations. The calibrated polarization of 3C286 of the three epochs were 8.9% ± 0.1%, 8.7% ± 0.2% and 9.0% ± 0.1% for polarization degrees and 30.4° ± 0.3°, 33.8° ± 0.5° and 29.4° ± 0.3° for polarization angles. Considering that the ionosphere can generate a Faraday rotation of 1°–3° in polarization angle at the L band33, our results were consistent with the intrinsic polarization degree of 9.5% and polarization angle of 33° of 3C286 at 1,450 MHz (ref. 34). In addition to the polarization observations of L1544 and 3C286, we observed the circularly polarized OH maser source IRAS02524+204635 to verify that our procedures produced consistent Blos, including the sign or direction of the magnetic field, as had been obtained previously.

The convolutions of the sidelobes of the Stokes V beam with the spatial gradient of the Stokes I emission may generate a false ‘S curve’ in the V spectrum27. To check the credibility of our Zeeman detections, we measured the Stokes V beam of FAST and convolved the beam with the Galactic Arecibo L-band Feed Array (GALFA) Stokes I cube36 of L1544. The convolved V spectrum showed a profile with a shape similar to the I spectrum and a strength less than 0.03% of the I spectrum, different from the ‘S curve’ patterns in the observed V spectrum. Meanwhile, the 19-beam receiver was rotated to −45°, 0° and 45° in the three epochs of the L1544 observations, and all of the three epochs showed ‘S curve’ patterns in the V spectra, indicating that our Zeeman results were true detections.

Although the data of the 19 beams of the FAST L-band receiver were simultaneously taken in our observations, only the polarization of the central beam was commissioned at the time of writing. The results represented in this work were made with only the central beam pointing towards the HINSA column density peak in Fig. 1. The Zeeman results of the 18 off-central beams will be published in the future.

Multiple Gaussians and radiative transfer fitting to I(v) and V(v)

We adopt the least-squares fits of multiple Gaussians with radiative transfer26 to decompose the I(v) into the HINSA, CNM and WNM components. The expected profile of I(v) consists of multiple CNM components providing opacity and also brightness temperature and a WNM component providing only brightness temperature:

$$I(v)={I}_{{\rm{CNM}}}(v)+{I}_{{\rm{WNM}}}(v).$$
(1)

The ICNM(v) is an assembly of N CNM components

$${I}_{{\rm{CNM}}}(v)=\mathop{\sum }\limits_{n=1}^{N}{I}_{{\rm{peak}},n}(1-{{\rm{e}}}^{-{\tau }_{n}(v)}){{\rm{e}}}^{-({\sum }_{m=0}^{M}{\tau }_{m}(v)+{\tau }_{0})},$$
(2)

where the subscript m with its associated optical depth profile τm(v) represents each of the M CNM clouds that lie in front of cloud n. The optical depth of the ith component is

$${\tau }_{i}(v)={\tau }_{i}{{\rm{e}}}^{-{[(v-{v}_{0,i})/{\sigma }_{v,i}]}^{2}}$$
(3)

in which τ0 represents the HINSA providing only opacity and no brightness temperature. For the WNM in the background

$${I}_{{\rm{WNM}}}(v)={I}_{{\rm{peak}},{\rm{WNM}}}{{\rm{e}}}^{-{[(v-{v}_{0,{\rm{WNM}}})/{\sigma }_{v,{\rm{WNM}}}]}^{2}}{{\rm{e}}}^{-{\sum }_{i=0}^{N}{\tau }_{i}(v)}.$$
(4)

The fitting of I(v) thus yields values for the intrinsic peak Stokes I emission (Ipeak), τ, v0 and the Gaussian dispersion (σv) of the components.

We consider the radiative transfer of V(v) in terms of right circular polarization (RCP) and left circular polarization (LCP). The Zeeman effect states that with the existence of Blos, the frequency of RCP shifts from its original frequency ν0 to ν0 + νz and the frequency of LCP shifts to ν0 − νz with νz = (Z/2) × Blos, where Z is the Zeeman splitting factor (2.8 Hz μG−1 for the H I 21-cm line). As the RCP and LCP are orthogonal components of radiation, the radiative transfer processes of RCP and LCP are independent of each other. For RCP, equation (1) becomes

$${T}_{{\rm{RCP}}}={T}_{{\rm{RCP}},{\rm{CNM}}}(v,{\tau }_{{\rm{RCP}},i})+{T}_{{\rm{RCP}},{\rm{WNM}}}(v,{\tau }_{{\rm{RCP}},i}),$$
(5)

where for the ith component, \({T}_{{\rm{RCP}},i}={I}_{{\rm{peak}},i}/2\), τRCP,i is the optical depth in the RCP radiation to substitute the τi in equation (3) with \({\tau }_{{\rm{RCP}},i}={\tau }_{i}({\nu }_{0}+{\nu }_{{\rm{z}},i})\) for Blos,i of the component, and the parameters of v0 and σv remain the same. Similarly, for LCP

$${T}_{{\rm{LCP}}}={T}_{{\rm{LCP}},{\rm{CNM}}}(v,{\tau }_{{\rm{LCP}},i})+{T}_{{\rm{LCP}},{\rm{WNM}}}(v,{\tau }_{{\rm{LCP}},i})$$
(6)

with \({T}_{{\rm{LCP}},i}={I}_{{\rm{peak}},i}/2\) and \({\tau }_{{\rm{LCP}},i}={\tau }_{i}({\nu }_{0}-{\nu }_{{\rm{z}},i})\). The fitting of \(V(v)={T}_{{\rm{RCP}}}-\)\({T}_{{\rm{LCP}}}+cI(v)\), which includes a c term accounting for leakage of I(v) into V(v), thus yields values for Blos of the components. In Extended Data Table 1, we list the parameters of the components obtained from least-squares fits to I(v) and V(v). The leakage of our HINSA Zeeman observations is c = 0.034%.

HINSA and H2 column density maps

L1544 is a low-mass prestellar core in the Taurus molecular cloud complex at a distance of about 140 pc. The core has a size of about 0.1 pc (ref. 37), presumably formed out of a parsec-long elongated molecular ridge7, which, for simplicity, we refer as the molecular envelope. We show the HINSA and H2 column density maps of L1544 in Fig. 1a, and we use the H2 column density map to calculate the NH2 of the envelope and core at the beams of FAST and Arecibo observations in Table 1. The HINSA column density map is a revision of Fig. 8 in ref. 4. To derive the H2 column density map, we retrieved the level-2.5 processed, archival Herschel images that were taken at 250 μm/350 μm/500 μm using the SPIRE instrument38 (observation ID 1342204842). We smoothed the Herschel images to a common angular resolution of the 36″ beam at 500 μm and regridded the images to the same pixel size of 6″. We performed least-squares fits of the 250 μm/350 μm/500 μm spectral energy distributions weighted by the squares of the measured noise levels to derive the pixel-to-pixel distributions of dust temperature Td and dust optical depth τν using \({S}_{\nu }={{\Omega }}_{m}{B}_{\mu }({T}_{{\rm{d}}})(1-{{\rm{e}}}^{-{\tau }_{\nu }})\), where Sν is the flux density at frequency ν, Ωm is the is the solid angle of the pixel, Bμ(Td) is the Planck function at Td and \({\tau }_{\nu }={\tau }_{230}{(\nu ({\rm{GHz}})/230)}^{\beta }\) with a dust opacity index β of 1.8. Next, we obtained the H2 column density with \({N}_{{\rm{H2}}}=g{\tau }_{230}/({\kappa }_{230}\,{\mu }_{{\rm{m}}}{m}_{{\rm{H}}})\), where g = 100 is the gas-to-dust mass ratio, κ230 = 0.09 cm2 g−2 (ref. 39) is the dust opacity at 230 GHz, μm = 2.8 is the mean molecular weight and mH is the atomic mass of hydrogen. To estimate the uncertainties in the H2 column density, we used a Monte Carlo technique. For each pixel, we created artificial 250 μm/350 μm/500 μm flux densities by adding the original flux densities with normal-distributed errors taking account the uncertainty in the measured flux and a 10% correlation for the calibration uncertainty in SPIRE40. We then estimated the uncertainty in each pixel with 1,000 fittings of the H2 column density. The NH2 and its uncertainty in Table 1 were obtained from the convolutions of the H2 column density map and uncertainty map with the FAST and Arecibo beams.

Note that the equivalent H2 column density NH2 of the CNM1 are derived from H I data towards 3C132 and 3C133, a method different to the NH2 of the L1544 envelope and core that are derived from dust emission. Therefore, in addition to the statistical errors listed in Table 1, there is a systematic difference between the NH2 derived from the two methods. Considering that the regime traced by dust emission can be different from those traced by HINSA or OH, which is particularly noticeable from the different spatial extents of dust and HINSA in Fig. 1a, we expect that the systematic difference could be as large as a factor of a few. As the values of λ between CNM1 and L1544 are different at least by a factor of 13, the systematic difference between the two methods should not change the qualitative conclusion of this work.

In Fig. 1a, the peak of HINSA column density appears to be shifted from the centre of L1544 by 0.15 pc and the 70% and 90% contours of the peak HINSA column density do not enclose L1544. We note that such offset has also been seen for other dense gas tracers in prestellar cores41. The core geometry may not be as simple as envisioned in idealized theories, where the dense core sits near the centre of a lower-density molecular envelope. In particular, the L1544 core appears to sit near one end of an elongated molecular (and dust) ridge, which roughly coincides with the region traced by the HINSA. Such an offset can result from complexities in chemistry and formation history, but does not affect the main science result of this work, namely, the HINSA Zeeman probes the magnetic fields of the current molecular ridge that is the progenitor of the dense core.

Maximum likelihood

We adopt the analysis of maximum likelihood18 to study the uniformity of magnetic fields in the envelope of L1544. Assuming that the true Blos follows a Gaussian distribution with mean B0 and intrinsic spread σ0, the likelihood lj for a single observation in a set of N measurements (j = 1, ..., N) to measure Bj with Gaussian error σj is proportional to the convolution of the probability \(\exp [-{(B-{B}_{0})}^{2}/2{\sigma }_{0}^{2}]/\sqrt{2{\rm{\pi }}}{\sigma }_{0}\) for the magnetic field to have a true value of B with the probability \(\exp [-{(B-{B}_{j})}^{2}/2{\sigma }_{j}^{2}]/\sqrt{2{\rm{\pi }}}{\sigma }_{j}\) of observing a value Bj of the field. Therefore, lj is the integral over all possible true values of the magnetic field

$${l}_{j}={\int }_{-\infty }^{\infty }{\rm{d}}B\frac{\exp [-{(B-{B}_{j})}^{2}/2{\sigma }_{j}^{2}]}{\sqrt{2{\rm{\pi }}}{\sigma }_{j}}\frac{\exp [-{(B-{B}_{0})}^{2}/2{\sigma }_{0}^{2}]}{\sqrt{2{\rm{\pi }}}{\sigma }_{0}}.$$
(7)

Although the overall likelihood \( {\mathcal L} \) for a set of observations is the product of individual likelihoods of the observations \(( {\mathcal L} ={\prod }_{J=1}^{N}{l}_{j})\), the B0 and σ0 can be estimated by maximizing the likelihood \( {\mathcal L} \). After performing the integration in equation (7) and some algebraic manipulations

$$ {\mathcal L} ({B}_{0},{\sigma }_{0})=\left(\mathop{\prod }\limits_{J=1}^{N}\frac{1}{\sqrt{{\sigma }_{0}^{2}+{\sigma }_{j}^{2}}}\right)\exp \left[-\frac{1}{2}\mathop{\sum }\limits_{j=1}^{N}\frac{{({B}_{j}-{B}_{0})}^{2}}{{\sigma }_{0}^{2}+{\sigma }_{j}^{2}}\right]$$
(8)

Extended Data Fig. 1 shows the distribution of \( {\mathcal L} \) as functions of B0 and σ0 and the probability distributions of B0 and σ0 by integrating \( {\mathcal L} \) along the B0 axis and the σ0 axis, respectively. The probability distribution of B0 is similar to a normal distribution with a mean value of +4.1 μG and a standard deviation of 1.6 μG. The probability distribution of σ0 is highly asymmetric as the values of σ0 cannot be negative. The first, second and third quartiles of the σ0 distribution are 0.6 μG, 1.2 μG and 2.4 μG. We therefore suspect that the Zeeman measurements in the L1544 envelope can be explained by a magnetic field with B0 = +4.1 ± 1.6 μG and \({s}_{0}={1.2}_{-0.6}^{+1.2}\,{\rm{\mu }}{\rm{G}}\).

Inclination angle of magnetic field

Given the uniformity of magnetic fields in the envelope of L1544 and CNM1 is well constrained by the maximum likelihood analysis, the coherent Blos suggests that the inclination angles of magnetic fields in the CNM1 and L1544 envelope are likely to be similar, or a special geometry of magnetic field structure across multi-scales and multi-phases of the interstellar medium is needed. In contrast, the Blos of the L1544 envelope and core differ by a factor of 2.6. There are two physical explanations for the 2.6-times difference between the HINSA and OH Zeeman measurements. First, the OH measurement probably samples a denser gas than the HINSA measurement, as the column density along the OH sightline is twice that along the HINSA sightline (Table 1). As the magnetic field strength in molecular clouds tends to increase with number density15, the stronger field is naturally expected in the denser core. Alternatively, the inclination angle of the L1544 core magnetic field could differ substantially from that of the coherent field. As we cannot rule out the second possibility, an assumption of similar inclination angles of the magnetic fields in the L1544 envelope and core thus is required to calculate the relative values of λ.

We note that dust polarization observation may give some clues as dust polarization traces the position angle of the plane-of-sky component of the magnetic field. The near-infrared polarization observations of L154442 indicate that the mean position angle of the magnetic field towards the core location of the Arecibo beam is 29.0°–36.9°, and the mean position angles of magnetic fields towards the four envelope locations of the GBT beams are 30.5°–55.8°. The difference in the position angles between the core and envelope thus may be about 10°–20°.

We perform Monte Carlo simulations43 to study whether the 2.6-times difference between the Blos of the L1544 envelope and core can be explained by different inclination angles. The simulations randomly generate two unit vectors in three dimensions, and then measure the difference between the inclination angles and the difference between the position angles of the two vectors. The probability of the cases that the line-of-sight length of one vector is 2.6-times larger than that of the other is roughly 0.19. For those cases, the mean difference between the inclination angles and the mean difference between the position angles of the two vectors are about 38° and 45°, respectively. As the probability of 0.19 is small and the difference of about 45° between the simulated position angles is about a factor of two to four times larger than the difference of about 10°–20° between the observed position angles, it is less likely that the 2.6-times difference between the Blos of the L1544 envelope and core can be solely explained by different inclination angles.

CCS Zeeman measurements

Ref. 20 reported a CCS Zeeman detection of 117 ± 21 μG in the dense core of TMC-1 that has an estimated H2 column density of 3 × 1022 cm−2, which is four-times higher than that probed by OH Zeeman measurements in L1544 and nearly one order of magnitude higher than that probed by our HINSA measurements. It appears to provide further support to the evolutionary scenario suggested by our HINSA measurements: namely, once the gas loses its magnetic support during the transition from the CNM to the molecular envelope (or ridge) and becomes magnetically supercritical, there is no longer any need to lose magnetic flux further (relative to the mass) for a piece of the envelope/ridge to condense into a (magnetically supercritical) core (for example, the L1544 core probed by OH) and for the core to evolve further by increasing its column density (for example, the TMC-1 core probed by CCS).

Technically, we note that one potential source of significant uncertainty in frequency shift, namely the uncertainty of beam squint, was not included in the CCS result, which may affect the level of significance. In comparison, the HINSA measurement is robust with a greater than 10σ significance with the beam squint and velocity gradient being taken into account by convolving the FAST Stokes V beam with the Stokes I cube of L1544 (see the third paragraph in the data reduction section in Methods).