MAGNETIC FIELD STRENGTH MAPS FOR MOLECULAR CLOUDS: A NEW METHOD BASED ON A POLARIZATION–INTENSITY GRADIENT RELATION

The Stokes I dust continuum contours and the magnetic field segments (red) from the dust polarization emission (Tang et al. 2009) are shown in Figure 1 for W51 e2. Interestingly, field segments tend to be perpendicular to Stokes I emission contours, or equivalently, polarization is preferentially tangential to the emission contour lines. Instead of further working with intensity contours, we alternatively choose to investigate the intensity gradients (blue thick segments in Figure 1). The top panel in Figure 2 displays the magnetic field position angles (P.A.s) versus the intensity gradient P.A.s. A close correlation becomes manifest. Both P.A.s are defined counterclockwise starting from the north in a range of 0°–180°. The few data points (black empty squares) which originally do not seem to follow the trend, belong to pairs where both segments are close to P.A. = 0 (vertical orientation), but with one of the P.A.s slightly rotated to the left and the other one slightly rotated to the right-hand side of the vertical. Despite being closely correlated in orientation (small absolute differences), these pairs would mimic a bad correlation simply due to the definition of the P.A.s. For these cases, the original range of the magnetic field P.A. (0°–180°) is extended by adding their relative difference in P.A.s to their intensity gradient P.A.s. For example, a pair with an intensity gradient and magnetic field P.A. of 175° and 5° has a relative difference of 10° and is therefore displayed as (175°,185°). This redefinition does not alter the correlation in any way, but is merely a consequence of the ambiguity in the conventional definition of the 180°P.A. range.

Zoom In Zoom Out Reset image size

We note that the correlation is not limited to those pairs directly in the core-collapsing region, but seemingly still holds for the most northern segments where the emission contours start to deviate from circular shapes. The correlation—in the definition of Pearson's linear correlation coefficient—is 0.95. It is important to remark that, despite this tight correlation, the distribution of the differences (between magnetic field P.A.s and intensity gradient P.A.s) is non-Gaussian with an absolute mean of ∼20°. A Kolmogorov–Smirnov test rejects the null hypothesis (differences following a Gaussian distribution) with a p-value of ∼10⁻⁷. As a matter of fact, the relative differences show a bimodal distribution (bottom panel in Figure 2), which results from a systematic clockwise and counterclockwise rotation of the magnetic field with respect to the intensity gradient. This further points toward a systematic origin of this correlation. All of this leads to the conclusion that magnetic field and intensity gradient cannot simply be radially aligned with the apparent misalignment just being the result of some measurement uncertainty. It is the focus of an ongoing study where this correlation is further analyzed on a larger data sample in the context of the evolution of molecular clouds (P. M. Koch et al. 2012, in preparation).

In the following section, this observed correlation and the measured non-Gaussian differences in position angles of the magnetic field and the intensity gradient are further used to investigate the field strength.

3.2.1. Local Field Strength

Molecular clouds are the sites of interactions of various forces. Gravity, magnetic fields, expanding shells such as from expanding H ii regions or supernovae, accretion, and outflows all add to the dynamics, and their relative importance can vary with the evolutionary stage of the star formation site. Generally, MHD equations—possibly for various molecules and dust with individual coupling terms—are needed to provide an accurate and complete description of the system. In the following, aiming at making use of dust polarization measurements, a single MHD equation for the dust component is adopted.⁶ Under the assumptions of negligible viscosity and infinite conductivity (ideal MHD case) we focus on the force equation:

$$\begin{equation} \rho \left(\frac{\partial }{\partial t}+\mathbf {v}\cdot \mathbf {\nabla }\right)\mathbf {v} = -\mathbf {\nabla }\left(P + \frac{B^2}{8\pi }\right)-\rho \nabla \phi +\frac{1}{4\pi } \left(\mathbf {B}\cdot \mathbf {\nabla }\right)\mathbf {B}, \end{equation} \tag{ 1 }$$

where ρ and $\mathbf {v}$ are the dust density and velocity, respectively. $\mathbf {B}$ is the magnetic field. P is the hydrostatic dust pressure. ϕ is the gravitational potential resulting from the total mass contained inside the star-forming region. $\mathbf {\nabla }$ denotes the gradient. As usual, the left-hand side in Equation (1) describes the resulting action based on the force terms on the right-hand side. These include the gradients of the hydrostatic pressure terms of the gas, the magnetic field, and the gravitational potential together with the magnetic field tension term (last term on the right-hand side).

In the following, an interpretation of Equation (1) is derived which can be matched to a polarization observation as in Figure 1. We start by noting that the magnetic field tension term, ${\sim} (\mathbf {B}\cdot \mathbf {\nabla })\mathbf {B}$, can be rewritten as

$$\begin{equation} \left(\mathbf {B}\cdot \mathbf {\nabla }\right)\mathbf {B} =B\frac{\partial B}{\partial s_B}\mathbf {e}_{s_B}+B^2\frac{\partial \mathbf {e}_{s_B}}{\partial s_B} =B\frac{\partial B}{\partial s_B}\mathbf {e}_{s_B}+B^2\frac{1}{R}\mathbf {n}, \end{equation} \tag{ 2 }$$

where a magnetic field line, $\mathbf {B}=B\mathbf {e}_{s_B}$, is directed along the unity vector $\mathbf {e}_{s_B}$ with the generalized coordinate s_B along the field force line. A change in direction, $\partial \mathbf {e}_{s_B} / \partial s_B$, measures the magnetic field curvature 1/R and is directed normal to $\mathbf {e}_{s_B}$ along the unity vector $\mathbf {n}$. In an analogous way, $(\mathbf {v}\cdot \mathbf {\nabla })\mathbf {v}$ can be reformulated with the velocity vector $\mathbf {v}=v\mathbf {e}_{s_v}$ directed along the unity vector $\mathbf {e}_{s_v}$.

Generally, the gradient of the isotropic field pressure, ∼∇(B²/8π), in Equation (1) is not necessarily directed along the field direction $\mathbf {e}_{s_B}$, but can also have a component orthogonal to it. Therefore, when combining Equations (2) and (1), a term ${\sim} B({\partial B}/{\partial s_{B_{\bot }}})\mathbf {e}_{s_{B_{\bot }}}$ remains which needs to be compared to the curvature term ${\sim} B^2({1}/{R})\mathbf {n}$. Both $\mathbf {n}$ and $\mathbf {e}_{s_{B_{\bot }}}$ are orthogonal to $\mathbf {e}_{s_B}$, but they are not necessarily collinear. Assuming $\partial s_{B_{\bot }}$ and R to be given approximately by the resolution of an observation, the significance of this additional term can be estimated as $| B({\partial B}/{\partial s_{B_{\bot }}})|/(B^2({1}/{R}))\approx {\Delta B_{\bot }}/{B}\ll 1$. Thus, the term can be omitted if the local change in the normal field component, ΔB_⊥, is small compared to the total field strength. This holds for any smooth and slowly varying function of field strength. Equation (1) can then be transformed into

$$\begin{equation} \rho v \frac{\partial v}{\partial s_v}\mathbf {e}_{s_v} +\rho v^2 \frac{\partial \mathbf {e}_{s_v}}{\partial s_v} = -\frac{\partial P}{\partial s_P}\mathbf {e}_{s_P} - \rho \frac{\partial \phi }{\partial s_{\phi }}\mathbf {e}_{s_{\phi }} +\frac{1}{4\pi }B^2\frac{1}{R}\mathbf {n}, \end{equation} \tag{ 3 }$$

where the resulting directions of the gradients of pressure and gravitational potential are given with $\mathbf {e}_{s_P}$ and $\mathbf {e}_{s_{\phi }}$, respectively. In the above formulation, generalized coordinates s_v, s_P, s_ϕ, s_B along the directions of the unity vectors $\mathbf {e}_{s_v}, \mathbf {e}_{s_P}, \mathbf {e}_{s_{\phi }}, \mathbf {e}_{s_B}$ are used. Generally, these vectors are not collinear, but have different three-dimensional directions with the condition that the vector sum of the right-hand side adds up to the resulting direction on the left-hand side. Each term in Equation (3) is intentionally written with its corresponding unity vector, defining its specific direction. This will later be used when interpreting Equation (3) together with an observed map. Equivalently, all unity vectors can be expressed in any regular coordinate system (e.g., spherical coordinates r, θ, ϕ) should this be of advantage. The unity vectors $\mathbf {e}_{s_v}, \mathbf {e}_{s_P}, \mathbf {e}_{s_{\phi }}, \mathbf {e}_{s_B}$ can be interpreted as volume-averaged directions, where the size of a volume element will depend on the resolution of an observation. In the above derivation we have neglected the partial time derivative, ∂/∂t, assuming stationarity.

Equation (3) describes the basic interaction of forces and the resulting motion. Linking some of these terms (in direction and strength) with observations will then allow us to solve for others. Our goal is to isolate the magnetic field strength B in the term $({1}/{4\pi })B^2({1}/{R})\mathbf {n}$. Therefore, the remaining terms in Equation (3) need to be identified and characterized. Dust in star formation sites is reacting to all the force terms in Equation (3). Maps of observed dust emission are reflecting the overall result of gravity, pressure, and magnetic field. Consequently, they represent a measure of the resulting motion, the left-hand side in Equation (3). We make the fundamental assumption that a change in emission intensity is the result of the transport of matter driven by the combination of the above-mentioned forces. Adopting this, it then follows that the gradient in emission intensity defines the resulting direction of motion on the left-hand side in Equation (3). For the directions of the gradients of the pressure and the gravitational potential we assume here, for simplicity, a spherically symmetrical molecular cloud, where the center can be identified from the peak emission. This assumption is relaxed in Section 4, where the method is generalized to an arbitrary cloud shape.

It is important to remark that the derivation so far is generally valid for three dimensions. Projection effects in both the integrated Stokes I and the polarized emission inevitably are present in an observation. In the following, we first proceed by identifying Equation (3) with an observed map in two dimensions. Projection effects are later addressed in Section 6.2. Figure 3 illustrates the further steps. The derivation of B is shown for a measured polarization (red) and intensity gradient direction (blue) in the first quadrant (with respect to the gravity center which here is supposed to coincide with the emission peak). For each measured pair of polarization and intensity gradient, its distance and direction from the gravity center defines the gravitational pull, $\rho ({\partial \phi }/{\partial s_{\phi }})\mathbf {e}_{s_{\phi }}$. Generally, gravity, magnetic field, and the intensity gradient directions are different. The schematic in Figure 3 illustrates how the vector sum of the various terms in Equation (3) is constructed. Solving for the magnetic field term then relies on measurable angles in the orientations between polarization and the intensity gradient (α) and the difference between the gravity and the intensity gradient directions (ψ). Applying a sinus theorem to the closed triangle, (∇P + ρ∇ϕ)/sin α = (1/4π)B²(1/R)/sin ψ, leads to the expression for the magnetic field strength:

$$\begin{equation} B=\sqrt{\frac{\sin \psi }{\sin \alpha }\left(\nabla P+\rho \nabla \phi \right)4\pi R}, \end{equation} \tag{ 4 }$$

where, in the most general case, all variables are functions of positions in a map.

Zoom In Zoom Out Reset image size

We note that—although illustrated for a simple, close to spherically symmetric case in Figure 3 with ∇P and ρ∇ϕ aligned and pointing toward the same center—Equation (4) is generally valid for any directions of ∇P and ρ∇ϕ. Thus, the above derivation is not restricted to a spherical collapsing core, but can be applied to any configuration where the various local force directions can be identified (see also Sections 4 and 6.1).

3.2.2. Local Field Significance

As an important further outcome of the method, the angle factor sin ψ/sin α in Equation (4) has a direct physical meaning. With the magnetic field tension force term F_B = B²/(4πR) and the gravitational and pressure forces |F_G + F_P| = |ρ∇ϕ + ∇P|, Equation (4) can be reformulated as

$$\begin{equation} \left(\frac{\sin \psi }{\sin \alpha }\right)_{{\rm local}}=\left(\frac{F_B}{|F_G+F_P|}\right)_{{\rm local}} \equiv \Sigma _B, \end{equation} \tag{ 5 }$$

where we have introduced Σ_B to define the field significance. This directly quantifies the local impact of the magnetic field force compared to the other forces involved (gravity and pressure gradient). Remarkably, the ratio is free of any model assumptions and simply relies on the two measured angles. It is purely based on the geometrical imprint by the various forces left in the field and intensity gradient morphologies. The angle factor, thus, also provides a model-independent qualitative criterion whether the magnetic field can prevent an area in a molecular cloud from gravitational collapse (Σ_B > 1) or not (Σ_B < 1). The field significance Σ_B and its implications for the mass-to-flux ratio and the star formation efficiency are investigated in Koch et al. (2012).

We briefly illustrate the expected features of our method when applied to a representative field and intensity configuration in a molecular cloud (Figure 4). A summary is given in Table 1. Field lines like (I) and (II) are expected in a core collapse when gravity is dragging inward the magnetic field which is flux-frozen to the dust particles. In the case of (I), the field line is close to the pole, little bent, and the gas is mostly moving along the field. Intensity gradients are mostly radial, pointing toward the center of gravity. Therefore, in the outer region (a), R ≈ const, δ small (sin α ≈ 1), ψ small (sin ψ ≈ ψ ≈ const), and thus Σ_B ≈ const. The field strength B will scale as ∼r^−κ, where κ > 0 describes a power-law profile resulting from density and gravitational potential gradient (Section 5.3). Small (random) changes in angles are possible due to local turbulence and can lead to some scatter. In the inner part (b), R might slightly decrease with the other parameters remaining roughly constant. This will lead to similar scalings for Σ_B and B. A field line of type (II) is expected around the accretion direction. The outer field (a) is still little bent (R ≈ const), δ ≈ const and small, ψ small, and therefore Σ_B ≈ const and B ∼ r^−κ as in (Ia). In (b) the field lines start to bend more significantly (R decreases) and δ increases (sin α decreases). ψ is small (≈ const) or possibly increasing if a non-spherical central core is being formed. Thus, Σ_B can increase as ∼1/sin α. The field strength is additionally modulated by $\sqrt{{R}/{\sin \alpha }}$, which can lead to a steeper or shallower profile than r^−κ depending on the exact values of R and α. In the most inner section (c), the field radius is minimized, field and intensity gradient directions are close to orthogonal to each other (δ ∼ π/2 modulo some turbulent dispersion, α ∼ 0), and ψ is small. The force ratio Σ_B can then possibly further increase locally. The field strength is thus mostly controlled by $\sqrt{{1}/{\sin \alpha }}$ and the power-law exponent κ, which will lead to highest values here. Overall, for (I) and (II) we thus expect the field strength to generally increase with smaller radius toward the center.

Zoom In Zoom Out Reset image size

Table 1. Qualitative Analysis

	ψ	α	R	ΣB	B	Comments
General
	↓↑	≠ 0	R	↓↑	↓↑, ${\sim} \sqrt{\sin \psi \cdot R}\cdot r^{-\kappa }$	ΣB and B↑ with larger angle ψ
	≠ 0	↓↑	R	↑↓	↑↓, ${\sim} \sqrt{\frac{R}{\sin \alpha }}\cdot r^{-\kappa }$	ΣB and B↑ with smaller angle α
Specific Cases
(I) (a)	small, ≈const	${\sim} \frac{\pi }{2}$	≈const	≈const	≈const· r−κ	Field line close to pole, ΣB ≈ const, B↑ with smaller r
(b)	small, ≈const	${\sim} \frac{\pi }{2}$	≈const	≈const	≈const· r−κ
(II) (a)	small, ≈const	${\sim} \frac{\pi }{2}$	≈const	≈const	≈const· r−κ	Field line around accretion direction, ΣB ≈ const, B↑ with smaller r
(b)	small, ≈const	↓	↓	${\sim} {\frac{1}{\sin \alpha }}$	${\sim} \sqrt{\frac{R}{\sin \alpha }}\cdot r^{-\kappa }$	ΣB↑, B↑ with smaller r steeper or shallower than r−κ (depending on R and α)
(c)	small, ≈const	↓, ∼0	↓, Rmin	${\sim} {\frac{1}{\sin \alpha }}$	${\sim} \sqrt{\frac{R}{\sin \alpha }}\cdot r^{-\kappa }$	ΣB↑, maximum field strength B
(III)	↑↓	≈const	≈const	↑↓	↑↓, ${\sim} \sqrt{\sin \psi }\cdot r^{-\kappa }$	Outer area, no clear signatures yet, B generally small due to r−κ
(IV)	large	≈const	≈const	large	${\sim} \sqrt{\sin \psi }\cdot r^{-\kappa }$	Disconnected area, no clear gravity center localized, B small due to r−κ
Limits
	→ 0	≠ 0	(i) Ra, (ii) → ∞b	→ 0	(i) → 0a, (ii) → ∞b	Triangle degenerates to single line
	≠ 0	→ 0	(i) R, (ii) → ∞	→ ∞	(i) → ∞, (ii) → ∞	Method fails, triangle not closed
	→ 0	→ 0	(i) R, (ii) → ∞	→1	(i) ${\sim} \sqrt{R}$, (ii) → ∞	Method fails, triangle not closed and degenerated to single line

Notes. Qualitative analysis of the force ratio Σ_B = sin ψ/sin α and the magnetic field strength $B=\sqrt{({\sin \psi }/{\sin \alpha })\rho \nabla \phi \cdot 4\pi R}$ (Section 3.3). The pressure gradient is neglected here. The profile of the density and the resulting gradient of the gravitational potential is represented with r^−κ with κ > 0 (Section 5.3). The specific cases labeled (I) to (IV) refer to Figure 4. The cases (Ib), (IIb), and (IIc) are not yet observed due to the lack of higher resolution data. Down arrows (↓) and up arrows (↑) indicate decreasing and increasing values, respectively. "≈const" indicates that small local variations can occur which lead to some scatter in Σ_B and B (Figure 11). Failure of the method—triangle cannot be closed—refers to Figure 3. ^aThe case ψ → 0 states generally that the intensity gradient does not deviate from the gravity direction for any orientation of the field lines. Therefore, the field strength must be insignificant, i.e., B → 0, which is also consistently reproduced by the method. ^bThe case α ≠ 0 ∼ π/2 with R → ∞—i.e., straight B field lines aligned with the intensity gradient—represents the strong field flux freezing situation (see Section 6.1) and consistently yields a dominating field strength (B → ∞, Σ_B → ∞).

Download table as:  ASCIITypeset image

Additional features are seen in Figure 1, which are possibly not directly linked to the above cases. Region (III)—corresponding to the NE patch around the offset (−0.3, 1.8)—shows a group of parallel B field segments with parallel intensity gradients. In this outer region the collapse is likely initiated, but does not yet reveal clear signatures. ψ is changing, δ ≈ const (sin α ≈ const), R ≈ const, and Σ_B can therefore increase or decrease. The field strength will generally be small here because of the low density. The relative field significance can still be important here, as gravity has not yet fully taken control to shape and align intensity gradient and field direction. Area (IV) corresponds to the most northern part in Figure 1. This region is possibly disconnected from the main core, likely belonging to another smaller core being formed here. The Σ_B values are dominated by large changes in ψ when linked to the gravitational center of the main core. Due to the lack of a clear identification of a local gravity center (in which case ψ would be much smaller), Σ_B is probably overestimated for the two most northern segments. Possible solutions to this issue are addressed in Section 4.2. The field strength is again small here due to the low density.

From the derivation in Section 3, both angles α and ψ, the gravitational potential and the pressure P are generally functions of positions (x, y) in a projected plane of sky map. The magnetic field line radius R (curvature C ≡ 1/R) will generally also change with positions. It is, however, not possible to define a curvature based on a single isolated polarization segment. For polarization observations with a sufficiently large number of segments (as, e.g., in Figure 1), a curvature resulting from two adjacent (or more) polarization segments can be calculated.

Two neighboring magnetic field segments, separated by distance d, are assumed to be tangential to a field line of radius R. The difference in their position angles, ΔP.A., is a measure of how much the orientation of the field line has changed over d. The curvature C can then be written as (Figure 5)

$$\begin{equation} C\equiv \frac{1}{R}=\frac{2}{d} \cos \left(\frac{1}{2}\left[\pi -\Delta {\rm P.A.}\right]\right). \end{equation} \tag{ 6 }$$

With this definition of a local curvature, the calculation of the field strength B in Equation (4) is no longer limited to an average field curvature, but can be precise to any locally varying field structure.

Zoom In Zoom Out Reset image size

We remark that the above definition for C is not necessarily unique. In particular, determining which segments connect to a field line can be non-trivial. Additionally, for some neighboring segments, as can also be seen in Figure 1, it is obvious that they cannot be both tangential to the same field line. The curvature as defined in Equation (6) is therefore rather a local-averaged curvature. In this view, a local curvature can not only be defined with two segments which are roughly lining up, but also with two segments side by side. The local field radius map in Figure 6 (middle left panel) is calculated by averaging the local curvatures between two adjacent segments over all available closest neighbors. A closest neighbor is defined through the grid spacing of the map. Depending on the resolution and the detection of a polarization structure, the definition of C can possibly be adjusted or refined.

Zoom In Zoom Out Reset image size

For W51 e2, the local field radius derived from Equation (6) is in the range between ∼1'' and ∼17'' (Figure 6, middle panels). Whereas most areas in the map show rather uniform patches with R ∼ 1''–3'', regions in the northeast and the southwest around the accretion plane show larger radii. This is possibly a consequence of the accretion flow which tends to align field lines until some limiting inner radius where the field will bend. Yet more likely, this is a boundary effect because only two almost parallel segments are available here for averaging. This leads to unphysically large radii. The local radius averaged over the map is 〈R〉 ≈ 33. Discarding the extreme values gives 〈R〉 ≈ 1''.

Unless the gas distribution in a molecular cloud is (close to) perfectly symmetrical in azimuth around an emission peak, ∇ϕ will differ in direction and strength from a simple radius-only dependence in the gravitational potential. Thus, for an arbitrary mass distribution, the direction and strength of the local gravitational pull, $\rho ({\partial \phi }/{\partial s_{\phi }})\mathbf {e}_{s_{\phi }}$ in Equation (3), can no longer simply be constructed with the direction and distance to the reference center (emission peak in Figure 1). Instead, at any position in a map, the distribution of all the surrounding mass has to be taken into account to calculate the resulting local gravitational force. Since the method in Section 3 is based on measuring angles between various force directions (including gravity) at all positions in a map, any substantial local deviations from a simple global gravity model can be relevant.

For a discrete mass distribution (m_i), the resulting total potential ϕ at a location $\mathbf {r}$ is: $\phi (\mathbf {r})=-G\sum ^n_{i=1} {m_i}/{|\mathbf {r}-\mathbf {r}_i|}$, where G is the gravitational constant. We assume that a detected dust emission (like in Figure 1) is proportional to the gravitating mass, and that it further traces reasonably well the overall gravitating potential. Additionally, it is assumed that the dust temperature is roughly constant over the emission map. In this case, there is no weighting for m_i when linking it to the integrated dust emission I_i. The local gravity force (per unit mass), $\mathbf {g}(\mathbf {r})=-\nabla \phi (\mathbf {r})$ is then

$$\begin{equation} \mathbf {g}(\mathbf {r})=G\sum ^n_{i=1}\frac{I_i\cdot f_i}{|\mathbf {r}-\mathbf {r}_i|^2} \cdot \mathbf {e}_i, \end{equation} \tag{ 7 }$$

where f_i is a conversion factor linking dust emission and total mass, and $\mathbf {e}_i$ is the direction between each detected position with dust emission ($\mathbf {r}_i$) and the local position ($\mathbf {r}$) where the gravity force is evaluated.

We remark that, if f_i ≡ f over the entire map, the local gravity direction can be derived without prior knowledge of f. For a typical observation, the mass (intensity) distribution is obtained by pixelizing a map with the beam resolution. Figure 6 (bottom left panel) shows the deviations from an azimuthally symmetrical potential for W51 e2 when calculating the local gravity directions following Equation (7) for f_i ≡ f. This analysis shows that, although W51 e2 (Figure 1) appears to be a rather spherically symmetrical collapsing core, the deviations are systematically positive or negative on one or the other side of the hourglass axis. Minimum deviations are roughly along the southeast–northwest direction, which is consistent with the hourglass axis. The average absolute deviation is ≈9°. As further discussed in the Sections 6.4 and 6.5, already this simple symmetry analysis points toward a systematic influence of the magnetic field.

These deviations provide a correction for the angle ψ (angle between the gravity direction and the intensity gradient, Figure 3), which originally is calculated assuming an azimuthally symmetrical potential (Figure 7, top left panel). Consequently, also the magnetic field strength B in Equation (4) needs to be recalculated. The angle α, as defined relative to the intensity gradient, remains unchanged with this local gravity correction.

Zoom In Zoom Out Reset image size

We finally note that further refinement of the local gravity correction introduced here is possible. Local overdensities—which follow the flow of the global dynamics—can eventually locally collapse once they pass the Jeans mass. The northwest extension in Figure 1 is possibly such a case. Using Equation (7), however, the resulting local gravity direction will always be dominated by the main mass (emission). Consequently, the angle ψ might be over- or underestimated. This defect can be corrected by identifying local collapse centers with a local maxima search and a suitable threshold criterion.

5.1.1. sin ψ-map

The angle ψ measures the difference between the orientations of the intensity gradient and gravity. Small values in ψ (or sin ψ) therefore indicate that a change in the local intensity structure is closely following global gravity. The top panels in Figure 7 show a map of sin ψ together with its trends as a function of radius and azimuth. The peak emission (Figure 1) is again assumed as the reference center, and the gravity direction is simply taken as the radial direction from the center to any measured intensity location. Small values of sin ψ—in the range between 0 and 0.2—are found for most of the magnetic field segments pointing toward the center. An X-pattern is revealed made up from the majority of the segments. There is a trend of larger sin ψ values in the northeast and the southwest where likely accretion is ongoing. The magnetic field segments in the far north show sin ψ values close to one. Here, the intensity gradient directions are very different (by up to 90°) from the global gravity directions. This northern patch, therefore, seems to be decoupled from the main collapsing core (Section 3.3, case (IV)). sin ψ averaged over the entire map gives 〈sin ψ〉 ≈ 0.29 (〈ψ〉 ≈ 20°). The average decreases to about 0.19 (〈ψ〉 ≈ 12°) if the six polarization segments in the far north are ignored. A hint of an increase in sin ψ with radius is apparent in the top right panel in Figure 7. Again, the most northern patch dominates values beyond a distance of ∼15. No clear structure in azimuth is revealed.

Since ψ measures a deviation from gravity, mass distributions differing from an azimuthal symmetry can introduce corrections to the local gravity direction (Section 4.2). The top left panel in Figure 8 displays the sin ψ-map calculated from the generalized gravitational potential in Equation (7). The overall patterns are still similar to the uncorrected sin ψ-map, but with smaller values in the accretion areas. As a result, the average is reduced to 〈sin ψ〉 ≈ 0.19 (〈ψ〉 ≈ 12°) and to about 0.09 (〈ψ〉 ≈ 6°) when discarding the most northern patch.

Zoom In Zoom Out Reset image size

5.1.2. sin α-map

The angle factor (sin ψ/sin α)_local = (F_B/|F_G + F_P|)_local ≡ Σ_B measures in a model-independent way the local impact of the magnetic field force (Section 3.2.2). The top panels in Figure 6 illustrate the results. Values for the relative field significance, Σ_B, range from ∼0.1 to ∼1.3, with an average over the entire map 〈Σ_B〉 ≈ 0.33. Thus, the magnetic field can balance only about one-third of the gravitational force. Consequently, it cannot prevent the core from any further collapse. We stress that this result is obtained without any assumption about the mass involved. The top panels in Figure 6 further show that the impact of the magnetic field changes as a function of position. Generally, the field influence seems to become weaker toward the center where gravity likely becomes more and more dominant. This result is further explored in Koch et al. (2012).

Maps for the magnetic field strength B are presented in this section for two scenarios: first, only the newly introduced force ratio Σ_B is taken into account, leaving the other parameters in Equation (4) constant. This will illustrate the influence of Σ_B. Second, the field strength is calculated by combining Σ_B together with an analysis of the density profile and the resulting gravitational potential. A first-order estimate of the field strength, comparing the gravitational force F_G = GM_Rρ/R²_G and the magnetic field tension F_B = B²/(4πR_B) is given in Tang et al. (2009). Since W51 e2 is a collapsing core, F_G > F_B, an upper limit of the field strength, B < 19 mG, can be derived. M_R(≈ 220 M_☉) refers to the gas mass enclosed within a radius R_G(≈ 1''). The radius of a magnetic flux tube is R_B(≈ 1''). The mass density ρ is estimated with the gas volume number density $n_{{\rm H}_2} (\approx 2.7 \times 10^7$ cm⁻³). G is the gravitational constant. The above value, B < 19 mG, will serve as a reference for the magnetic field strength maps presented in the following. We remark that Equation (4) is consistent with the above estimate, $B=\sqrt{4\pi R_B G ({M_R \rho }/{R_G^2}})=\sqrt{4\pi R_B \rho \nabla \phi }$, when setting sin ψ/sin α = 1 and neglecting the pressure gradient ∇P.

In the first scenario, the B map for W51 e2 is calculated by taking into account only the angles ψ and α. Field radius, potential, and density are set constant as above. Assuming local changes in temperature and density to be negligible compared to gravity, the pressure gradient ∇P is omitted from Equation (4) in the following. From the measured ranges in ψ and α (Figure 7; sin ψ ≈ 0.1–0.9, sin α ≈ 0.7–1), it is then obvious that mostly ψ will dominate the range and structure in a magnetic field strength map (Figure 9). Values of the factor sin ψ/sin α = Σ_B range from ∼0.1 to ∼1.3, with an average over the entire map 〈Σ_B〉 ≈ 0.33 (Figure 6, top panels). The resulting patterns are very similar to the ones in the sin ψ-map (Figure 7, top left panel). Resulting variations in B are then from a few mG to ≈20 mG. The map average value is 〈B〉 ≈ 10 mG. Compared to the reference value above, the correction factor sin ψ/sin α generally leads to a lower field strength. Discarding the most northern patch—where B is likely overestimated due to the lack of a clear local gravity center (Section 3.3)—the value becomes 〈B〉 ≈ 8 mG. Obviously, when leaving curvature, gravitational potential, and density constant, the radial profile of the field strength is determined by Σ^1/2_B. As the field-to-gravity force ratio decreases toward the center, consequently the field strength also decreases and its largest values tend to be in the outer zones (Figure 9). This result needs to be compared with the following complete analysis.

Zoom In Zoom Out Reset image size

In the second scenario, the field strength is calculated by additionally taking into account the profiles of density and gravitational potential. As discussed later in Section 6.2, the force ratio Σ_B is little or not at all affected by projection effects. Therefore, we proceed to reconstruct three-dimensional radial profiles assuming spherical symmetry. For that purpose, the azimuthally averaged dust emission profile from Figure 1 is deprojected with an Abel integral. The resulting three-dimensional radial profile (r) is well described with a functional form f(r) ∼ 1/(r + 0.35)³, whereas the initial projected map (R) is well fit with a profile ∼1/(R + 0.25)². The fitting is performed on profiles with values binned to half of the synthesized beam resolution. Assuming the total mass of the e2 core to be distributed similar to its dust component, the cumulative mass m(r) ∼ ∫^r₀f(r') 4πr'² dr' and the resulting gradient of the potential ∼m(r)/r² can be derived. Combining all together (Equation(4)), the magnetic field strength profile scales as B(r) ∼ (Σ_B(r) · f(r) · m(r)/r²)^1/2. Whereas the dust emission map allows us to derive this functional form, its absolute normalization remains to be determined. Adopting 〈Σ_B〉 ≈ 0.33 (Section 5.2)—i.e., the observed averaged magnetic field strength can balance only about one-third of the gravitating mass M_R(≈ 220 M_☉)—we normalize B so that 〈B²〉 = 〈Σ_B〉 · (19 mG)², where the latter field strength is equivalent to 220 M_☉ (Tang et al. 2009). This gives 〈B〉 ≈ 9.4 mG. The final field strength map and its radial and azimuthal trends are shown in the left panels in the Figures 10 and 11. An increase in field strength with smaller radius becomes apparent. A profile ∼r^−1/2 is overlaid for illustration. A hint of an increase in field strength around r ≈ 15 possibly points toward the forming core in the northwest. This is found by still assuming a dominating gravity center at the emission peak. It is obvious that the first map (Figure 9) cannot represent the correct field strength, but merely reflects the relative field significance. Only by adding profiles for the density and the gravitational potential can the magnetic field profile correctly be derived. Density and gravitational potential only (setting Σ_B ≡ 1) would lead to a field profile ∼r^−1.7. Taking into account the decreasing field-to-gravity force ratio with radius, leads to the shallower scaling B(r) ∼ r^−1/2. With the density profile, ρ ∼ r⁻³, this connects field and density as B(r) ∼ ρ^1/6. Whereas the dust emission shows little variations in azimuth (close to spherical symmetry), the local position-dependent measurement of Σ_B shows some variations. Consequently, features in the magnetic field strength map start to be revealed, with possibly regions of increased strength in the accretion zones in the northeast and the southwest directions. In order to trace the magnetic field to the very center, higher resolution observations will be needed.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We remark that the projected dust emission can also reasonably well, although with larger residuals, be fit with a ∼R⁻¹ profile. The deprojected profile f(r) ∼ 1/(r + 0.25)² then gives a magnetic field profile B(r) ∼ r^−0.3 and B(r) ∼ ρ^3/20 (not shown). This leads to about 15% lower values in the field strength in the center with a similar mean field strength 〈B〉 ≈ 9.8 mG. All the other features are left unchanged. In the above final field strength map we have neglected the variations in the field radius and assumed it to be constant (Figure 6, top panels). This ignores two zones where large field radii are likely a consequence of boundary effects (Section 4.1).

Finally, the influence of the local gravity correction is illustrated in the right panels in the Figures 10 and 11. Generally, features tend to be smoother. The field strengths in the accretion zones are less pronounced and a single center becomes apparent. The average strength is about 7.7 mG. The radial profile of the binned values remains close to B(r) ∼ r^−1/2 with maximum values around 19 mG. Structures in azimuth remain similar although with smaller variations.

We stress that, although we are making use of the total gravitating mass M_R—which sets an upper limit for the field strength and defines its order of magnitude—only the newly derived force ratio Σ_B leads to a proper normalization and opens the possibility to derive a radial profile. This has become possible by solving the MHD force equation by the method introduced in Section 3.2.

Errors ΔB in the final B field maps are estimated by propagating typical measurement uncertainties through Equation (4). The individual polarization P.A. uncertainty of 5°–10° is assumed for Δα. A slightly smaller uncertainty, ∼3°–5°, can be expected for Δψ because of the averaging process in the interpolation when calculating the intensity gradient. This leads to ΔB ≈ ±(1/2B)(4πRρ∇ϕ) × 0.05, where the numerical factor follows from the above uncertainties combined with map-averaged values for 〈sin α〉 and 〈sin ψ〉. For any of the above scenarios and 〈B〉 ≈ 9 mG, we derive an average systematic error 〈ΔB〉_sys ≈ ±1 mG. Errors in the profiles for density and gravitational potential are small after binning and azimuthal averaging. They are, therefore, neglected here. Besides this systematic error, there is an additional statistical uncertainty in the field orientation resulting from the turbulent field dispersion. Adopting 5°–10° for this (Koch et al. 2010), it contributes to the final error budget with 〈ΔB〉_stat ≈ ±1 mG. We note that all uncertainties in angles are smaller than the mean value (≈20°) of the differences between magnetic field and intensity gradient P.A.s (Section 3.1). This further confirms that the correlation in Figure 2 cannot only result from systematic and statistical uncertainties.

In conclusion, density and gravitational potential profiles determine the order of magnitude of the field strength. The local features result from the angle factor sin ψ/sin α. These features remain robust with typical total errors (systematic and statistical) of ΔB ≈ ±(1–2) mG. The derived field strengths are also consistent with the level of magnetic fields typically detected by Zeeman effects on OH masers in compact H ii regions (a few mG up to ∼20 mG; Fish & Reid 2007).

The method presented in Section 3 leads to a local magnetic field strength. It, nevertheless, uses strictly speaking a not exactly local but larger scale property of the field, namely, its curvature (Figure 5). The shape of a field line is the result of various forces acting on the dust particles coupled to the magnetic field. These can be both long-ranging forces (like gravity) as well as locally acting perturbations (like turbulence). In order to solve Equation (1), we need to determine the larger scale curvature. The locally turbulent field can here act as a contaminant. This is in contrast to, e.g., the CF method which precisely uses the field dispersion around an average field direction, in combination with some turbulent velocity information, to derive a magnetic field strength. In principle, a curvature can be calculated from a minimum of two segments (Section 4.1), but the local turbulent field can introduce some uncertainty here. Assuming an isotropic turbulent field component, a local curvature can be derived by averaging over several neighboring curvatures. As long as reasonably smooth and connected polarization patterns are detected, this will reveal the larger scale field curvature. This is the case for the field morphology in Figure 1 where a clear signature of dragged-in field lines is apparent. Higher resolution data are then needed to trace the field in the inner most part of the collapsing cloud where the curvature is expected to change faster within a small area.

The method further relies on identifying a noticeable effect from the various forces in Equation (1) which shape the field on large scales. When assuming flux freezing, this means that B²/(8π) < P_gas + ρϕ + ρv². In this weak field case of flux freezing, the field lines are forced to move along with the gas (at least up to some scale where ambipolar diffusion sets in). As a consequence, gravity, intensity gradient, and field orientations can be individually determined which is needed to calculate the angles in Equation (4). In the strong field case, B²/(8π) ≫ P_gas + ρϕ + ρv², where matter can mostly move only along field lines, the method will yield ψ ∼ 0, α ∼ π/2, R → ∞, and therefore B → ∞. (See also Section 3.3.) Unless the pressure gradient ∇P can be identified, which then allows us to measure ψ ≠ 0, the method cannot further constrain B (see also Table 1). Geometrically, the method is based on closing a force triangle (Figure 3). Thus, if ψ ∼ 0 or α ∼ 0 or if both angles are zero at the same time, the triangle is no longer defined and the method fails.

In summary, whereas the CF method is based on the field dispersion but has to isolate a mean large-scale field direction, the method discussed here uses the large-scale field morphology shaped by various forces. Therefore, the method is more directly applicable to the weak field case. A successful extension to the strong field case is possible if pressure and/or gravitational potential gradients can be localized. Generally, coherent patches with smooth changes in field and intensity gradient directions are of advantage as this leads to a clear identification of the field curvature and the acting forces. The CF method needs both dust polarization continuum and molecular line velocity observations. The method here uses the dust continuum Stokes I and its polarized emission. Whereas the CF method provides one single statistically averaged value, the approach here leads to a map of field strength values.

Finally, we conclude this section by noting that the method here is not limited to dust continuum polarization. As long as a force triangle can be identified, the method is applicable to any measurement which reveals the large-scale field morphology in combination with gravity and/or pressure forces. Notably, molecular cloud polarization maps obtained from the GK effect, i.e., from spectral line linear polarization (Goldreich & Kylafis 1981, 1982), can be analyzed in the same way. This might open a way to study outflows which are often not observed in dust continuum due to their relatively small masses. The only complication possibly arises from the difficulty that it can be non-trivial to conclude whether the magnetic field is parallel or perpendicular to the polarization detected from the GK effect. In an even broader context, the method introduced here can also be applied to field morphologies detected from polarization by synchrotron emission or from Faraday rotation measures. The method is, therefore, not limited to molecular clouds, but can potentially also be adopted for, e.g., galaxy rotation measures.

The results in Section 5 were derived assuming the plane of sky projected geometry in Figure 3. Here, we address the question of how an original three-dimensional configuration affects the parameters in Equation (4) to calculate the magnetic field strength. The schematic in Figure 12 illustrates the connection. We note, that like in the two-dimensional case, gravitational pull (and pressure gradients), magnetic field tension, and the resulting intensity gradient have to form a closed triangle. The projected two-dimensional triangle can generally be deprojected with two different inclination angles: ι_a and ι_b, corresponding to the inclinations of the deprojected orientations of the magnetic field tension and the gravitational pull, respectively.

Zoom In Zoom Out Reset image size

Strictly speaking, projection effects are likely to be present to some extent in all variables in Equation (4). Having deprojected the intensity emission profile and adopted an average field radius (Section 5.3), we limit the discussion to the remaining factor sin ψ/sin α. In the following, lower indices "3" or "2" refer to three-dimensional or two-dimensional projected quantities. The deprojection of the side b₂ of the triangle (Figure 12) can be written as $b_3=\sqrt{\vphantom{A^A}\smash{\hbox{${b_2^2+\pi _b^2}$}}}$ and tan ι_b = π_b/b₂, where π_b is the length of the projection along the line of sight (z-direction). Analogous relations hold for the side a₂. Linking a sinus theorem in the three-dimensional triangle, sin ψ₃/sin α₃ = a₃/b₃, to the previously used two-dimensional expression sin ψ₂/sin α₂ = a₂/b₂, then yields

$$\begin{equation} \frac{\sin \psi _3}{\sin \alpha _3}=\frac{\sin \psi _2}{\sin \alpha _2} \cdot \sqrt{\frac{1+\tan ^2 \iota _a}{1+\tan ^2 \iota _b}} =\frac{\sin \psi _2}{\sin \alpha _2}\cdot \frac{\cos \iota _b}{\cos \iota _a}, \end{equation} \tag{ 8 }$$

where ι_a, ι_b ∈ [ − π/2, π/2]. It becomes obvious that for any inclination angles, the projection effects cancel out as long as gravitational pull (ι_b) and field tension (ι_a) are inclined by the same angles, i.e., ι_a ≡ ι_b. Furthermore, if ι_a and ι_b are different, projection effects are negligible as long as ι_a and ι_b are small. For ι_a and ι_b within about ±40°, the projection effect is within about ±20%. The left panel in Figure 13 displays the relevant function f(ι_a, ι_b) = cos ι_b/cos ι_a, which shows a large relatively flat saddle around one in the center. For extreme values, i.e., inclination angles |ι_a|, |ι_b| > 75°, the projection effect becomes severe (f → 0 or f → ∞). Besides having two totally different inclination angles ι_a and ι_b, we also investigate the projection effect as a function of the difference between the two inclination directions, Δ = ι_b − ι_a (Figure 13, right panel). Since the gravitational pull and the field tension directions are dynamically coupled, their inclination angles are likely to be not completely independent. The right panel in Figure 13 shows that projection effects are within about ±30% for a range in the differences Δ up to about 100° or more, up to inclinations of about ±45°. For larger inclinations f grows again. These situations, however, might not be of a concern, because at such large inclinations of the field lines, polarization emission (perpendicular to the field lines) might not be detectable. The analysis might, therefore, not be applicable.

Zoom In Zoom Out Reset image size

In summary, for identical inclination angles projection effects cancel out. The effect is small (up to ∼30%) or even negligible for one inclination up to ∼45° and the other one varying within about +50° and −50° (range in Δ up to about 100°). Only for extreme inclinations larger than about 75° projection effects can severely affect the method. Such situations might, however, not be detectable. It thus has to be emphasized that it is unlikely that severe projection effects are present if regular polarization patterns are observed. Consequently, the method introduced here leads to a total magnetic field strength.

Without any further modeling, the method does not include components such as, e.g., rotation or temperature variations. Effects due to rotation are not explicitly considered in Equation (1). Nevertheless, an additional term due to rotation on the right-hand side of Equation (1) would primarily lead to a change in the resulting direction of motion, $\rho v ({\partial v}/{\partial s_v})\mathbf {e}_{s_v} +\rho v^2 ({\partial \mathbf {e}_{s_v}}/{\partial s_v})$ in Equation (3), which was identified with the intensity gradient direction. Technically, this is equivalent to an additional inclination in the intensity gradient and can thus be approximated as a projection effect. Similarly, the orientations of the magnetic field lines dragged along with the rotating gas would be subject to an additional projection effect. As discussed in Section 6.2, the resulting uncertainty in the field strength is relatively small for a large parameter space of inclination angles. Alternatively, if the rotation velocity and axis are known, Equation (1) could be written in a rotating frame with an additional contribution from a centrifugal force.

Density and temperature combined make up for the detected dust emission. Significant variations in temperature could mimic mass concentrations which might then lead to wrong gravity directions and, therefore, to wrong magnetic field strength values. If concentrated locally such as, e.g., in the case of radiation pressure from a star, this is unlikely to alter the global features in the field strength over an entire map. Similarly, small shifts in the position of a dominating emission peak (like in Figure 1) introduce only minor variations in the field strength. Only if there are very irregular and significant changes in temperature over a map, is the method likely to fail. Additional information (such as from spectra) might then help to isolate the temperature effect.

Following up on Section 6.1, we expect the strong and weak field case of flux freezing to be present with different characteristics. In the strong field case—gas moving along the field lines—the intensity gradient direction is presumably closely aligned with the magnetic field orientation. On the other hand, in the weak field case—field lines dragged along and shaped by the forces acting on the gas—intensity gradient and field orientations can be at any angle. In the extreme case they can be at a 90° angle, when, e.g., a gas flow is bending a field line (case (IIc) in Figure 4). The angle δ(= π/2 − α), as previously introduced in Figure 3, measures the absolute difference in orientations between the intensity gradient and the magnetic field. Since the cases of weak and strong⁷ flux freezing categorize the relative importance of the field, we propose here the angle δ as a diagnostic measure for the dynamical role of the B field. Small angles (∼0) will point toward a dominating field, whereas larger and larger angles mean that the field is being overwhelmed by other forces.

The map in Figure 14 displays δ, revealing some areas with δ ≈ 0. These areas roughly coincide with the expected hourglass at larger radii. Here, the gas is possibly collapsing along the field lines. Additionally, it is likely that gravity has initially led to this field and intensity gradient configuration, with δ ≈ 0 being a mere consequence of that. Both at smaller radii toward the center and away from the hourglass toward the accretion plane, values for δ increase. This might be an indication that the field is dynamically less important here and that gravity has largely taken control. It needs to be stressed that this does not necessarily mean that the field strength is smaller here, because B further also depends on R, ∇ϕ, and ∇P (see also previous footnote). Moreover, the field can still provide pressure support which can slow down the gravitational collapse. This picture—with generally smaller δ values around the core and larger values toward the northwest corner—is consistent with the angle factor sin ψ/sin α (Section 3.2.2 and top panels in Figure 6), which identifies the field force to be less significant in the core. Figure 14 is necessarily complementary to the middle panels in Figure 7, because δ = π/2 − α.

Zoom In Zoom Out Reset image size

We have to note that the angle δ (like α) is subject to a projection correction. In the notation of Section 6.2: δ₃ = π/2 − α₃ with sin α₃ = sin ψ₃/sin ψ₂ · sin α₂ · f⁻¹(ι_a, ι_b). In the particular case of a close to spherical collapse, sin ψ₃/sin ψ₂ ≈ 1, because both the three-dimensional and the projected intensity show a clear azimuthal symmetry (Figure 1). Assuming f(ι_a, ι_b) to be approximately constant over a map (Section 6.2), then gives sin α₃ ≈ const · sin α₂. Thus, despite the projection effect, the dynamical importance of the field can still be understood from the relative differences in δ. The significance of the angle δ as a diagnostic measure is further analyzed on a larger data set together with the correlation from Section 3.1 in P. M. Koch et al. (2012, in preparation).

The method developed in Section 3 is based on identifying the three orientations of gravity, the magnetic field, and the emission intensity gradient. In particular, two angles are relevant in order to isolate and calculate the magnetic field strength B: the angle ψ between the gravity direction and the intensity gradient, and the angle α between the polarization orientation and the intensity gradient. Given the three orientations, an additional angle, ω, can be read off from a map. ω measures the difference in orientations between gravity and the magnetic field (Figure 3). Thus, it is a measure for how effectively gravity is able to shape the field lines and align them with the direction of the gravitational pull. Small angles (small values of sin ω) indicate close alignment. The bottom panel in Figure 7 shows sin ω for the case of W51 e2, assuming the center of gravity to be at the emission peak. The bottom panels in Figure 8 display again sin ω when a local gravity direction (Section 4.2) is adopted. The two figures are qualitatively similar. Applying a local gravity direction seems to sharpen the contrasts. The closest alignment (sin ω ≈ 0.1–0.2) is found at the outer radii close to the accretion plane in the west and in the northwest direction around the axis of the hourglass-like structure. An obvious poor alignment is in the northern area. Averaged over the entire map and excluding the northern patch we derive 〈sin ω〉 ≈ 0.21 (〈ω〉 ≈ 13°) and 〈sin ω〉 ≈ 0.14 (〈ω〉 ≈ 8°), respectively. Interestingly, after the local gravity correction, there is a very close resemblance between the sin ω-map and the sin δ-map (Figure 14). Since ω measures the magnetic field–gravity alignment and δ measures the magnetic field–intensity gradient alignment, the very similar maps, of course, indicate that the intensity structure is mostly resulting from gravity. This is also shown with the small values in the sin ψ-map (Figure 8, top panels). This is, of course, a consequence of the geometry in Figure 3 because all of the three angles are dependent: ω = π/2 − α − ψ with α = π/2 − δ. Nevertheless, the clear systematic differences in structures between the sin ω and the sin ψ-map (Figure 8) lead to the conclusion that the B field—at least in some regions—has kept its own dynamics and is not yet aligned and overwhelmed by gravity. This result is then also consistent with the finding of zones of weak and strong flux freezing (Section 6.4) which characterizes the dynamical role of the magnetic field. The systematic influence of the field is further manifest in the systematic deviations from an azimuthally symmetric emission structure, seen in the bottom panels in Figure 6. Projection effects (Section 6.2) are again likely to change only the absolute numbers, leaving relevant features and comparisons between maps still intact.