© ESO, 2011

1. Introduction

1.1. Motivation

Intergalactic magnetic fields are of great scientific interest. They are part of the intergalactic metabolism, they affect the propagation of heat and cosmic rays, and they are tracers of the dynamical state of the intergalactic medium. One way to retrieve their properties is to measure the Faraday rotation effect they cause on traversing polarised light. This has lead to significant insight into the presence, strength, and morphology of magnetic fields in galaxy clusters.

Although it is not possible to retrieve the exact magnetic field in position space from Faraday rotation data, the magnetic power spectrum is accessible and contains information on the typical field strength, correlation length, and the distribution of magnetic energy onto different length-scales. All these quantities can be predicted in magneto-genesis theories. The latter, the spectral energy distribution, is a sort of fingerprint of these theories. Magneto-genesis-scenarios in which the energy is injected on large scales by gas stirring typically exhibit Kolmogorov-like power law spectra (Subramanian et al. 2006; Enßlin & Vogt 2006). Whereas scenarios in which the magnetic fields are built up from small to large scales typically exhibit different spectral slopes (Schekochihin & Cowley 2006). Magnetic field measurements by Faraday rotation therefore have the potential to distinguish between these theories.

1.2. Faraday rotation

The plane of polarisation of light traveling to us through a magnetised plasma is changed by the wavelength-dependent Faraday rotation effect. If the Faraday rotating plasma is entirely located between the radio source (at physical depth z_s and at position x _⊥  in the plane of sky) and the observer (assumed to be at z = ∞), the effect is a plain change of the original polarisation angle χ₀ according to $\mathit{\chi }\mathrm{(}{x}\mathrm{\perp }\mathit{,\lambda }\mathrm{)}\mathrm{=}{\mathit{\chi }}_{\mathrm{0}}\mathrm{+}{\mathit{\lambda }}^{\mathrm{2}}\hspace{0.17em}\begin{array}{}\underset{\mathrm{􏽼}0\mathrm{􏽻􏽺}0\mathrm{􏽽}}{{\mathit{a}}_{\mathrm{0}}{\mathrm{\int }}_{{\mathit{z}}_{\mathrm{s}}}^{\mathrm{\infty }}\mathrm{d}\mathit{z}\hspace{0.17em}{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\mathrm{\perp }\mathit{,z}\mathrm{\right)}\hspace{0.17em}{\mathit{B}}_{\mathit{z}}\mathrm{\left(}{x}\mathrm{\perp }\mathit{,z}\mathrm{\right)}}\\ \mathrm{=}\mathit{\phi }\mathrm{(}{x}\mathrm{\perp }\mathrm{)}\end{array}\mathit{,}$(1)with a proportionality constant ${\mathit{a}}_{\mathrm{0}}\mathrm{=}{\mathrm{e}}^{\mathrm{3}}\mathit{/}\mathrm{(}\mathrm{2}\mathit{\pi }\hspace{0.17em}{\mathit{m}}_{\mathrm{e}}^{\mathrm{2}}{\mathit{c}}^{\mathrm{4}}\mathrm{)}$, with c, e, and m_e being the speed of light, the electron charge, and the mass, respectively. Furthermore, n_e(x) and B(x) are the electron density and the magnetic field within the Faraday rotating foreground medium, and φ(x _⊥ ) denotes the Faraday depth of the source at x _⊥ . Although the original polarisation angle χ₀ is unknown, measurements of the polarisation angle at a range of wavelengths λ provide this observationally, and also an estimator of the Faraday depth φ via a fit to Eq. (1). The observational estimate of the Faraday depth is called the rotation measure (RM), because it expresses the rate at which the linear polarisation rotates as a function of λ².

Some of the many studies of magnetic fields via Faraday rotation are Dennison (1979), Lawler & Dennison (1982), Dreher et al. (1987), Laing (1988), Hennessy et al. (1989), Taylor et al. (1990), Goldshmidt & Rephaeli (1993), Kim et al. (1991), Ge & Owen (1994), Feretti et al. (1995), Johnson et al. (1995), Johnson et al. (1995), Feretti et al. (1999a,b), Clarke et al. (2001), Taylor et al. (2001), Carilli & Taylor (2002), Taylor et al. (2002), Eilek & Owen (2002), Murgia et al. (2004), Laing et al. (2006), Xu et al. (2006), Taylor et al. (2007), Cho & Ryu (2009), Feain et al. (2009), and Battaglia et al. (2009), where we omitted works directly addressed in the following.

1.3. REALMAF

In this work, we present the code REALMAF (REAL MAgnetic Fields), which is aimed at measuring magnetic power spectra in galaxy clusters from maps of Faraday rotation measurements and apply this code to the Hydra A cluster.

As will be detailed in the following sections, our method relies on these informations and assumptions:

• the linearity of the measurement process (seeSect. 2.1);

• the maps of rotation measures and their errors (see Sect. 4.1);

• knowledge of the geometry of the Faraday active volume, which also depends on the projection angle θ of the jet and the cavity size (see Sect. 4.2);

• assuming that there is no polarised radio emission inside the Faraday active volume in front of the radio bubble (see Sect. 4.2);

• knowledge of the electron density profile n_e of the cluster atmosphere (see Sect. 4.2);

• assuming that the magnetic field strength B follows the electron density profile like $\mathit{B}\mathrm{\propto }{\mathit{n}}_{\mathrm{e}}^{\mathit{\alpha }}$ (e.g. Dolag et al. 2001,and see Sect. 3.4);

• assuming Gaussian probability distributions for the magnetic fields and the measurement errors in the rotation measure map (see Sect. 2.2);

• assuming statistical isotropy of B (see Sect. 2.3);

• assuming solenoidal fields (see Sect. 2.3);

• using a Bayesian approach to inference with a scale-independent prior for the power spectrum (see Sect. 3.1).

We use H₀ = 70 km s^-1 Mpc^-1, Ω_m = 0.3, Ω_Λ = 0.7 and a redshift of 0.0538 for Hydra A (de Vaucouleurs et al. 1991), which gives a linear to angular conversion factor 1′′ = 1.047 kpc. To follow the tradition in cosmic magnetism studies, we use cgs-units in this paper with the energy unit 1 erg = 10^-7 J and magnetic field strength unit 1 G = 10^-4 T.

2. Magnetic field statistics in Faraday screens

2.1. Measurement equation

The rotation measure differs from the Faraday depth as defined by Eq. (1) only by the measurement noise n: $\begin{array}{ccc}\mathrm{RM}\mathrm{\left(}{x}\mathrm{\perp }\mathrm{\right)}& \mathrm{=}& \mathit{\phi }\mathrm{\left(}{x}\mathrm{\perp }\mathrm{\right)}\mathrm{+}\mathit{n}\mathrm{\left(}{x}\mathrm{\perp }\mathrm{\right)}\\ & \mathrm{=}& \sum _{\mathit{i}\mathrm{\in }\mathrm{\left\{}\mathit{x,y,z}\mathrm{\right\}}}\hspace{0.17em}{\mathrm{\int }}_{{\mathit{z}}_{\mathrm{s}}}^{\mathrm{\infty }}\begin{array}{}\underset{\mathrm{􏽼}0\mathrm{􏽻􏽺}0\mathrm{􏽽}}{\mathrm{d}\mathit{z}\hspace{0.17em}{\mathit{a}}_{\mathrm{0}}\hspace{0.17em}{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\mathrm{\perp }\mathit{,z}\mathrm{\right)}\hspace{0.17em}{\mathit{\delta }}_{\mathit{i}\hspace{0.17em}\mathit{z}}}\\ {\mathit{R}}_{{x}\mathrm{\perp }\hspace{0.17em}\mathrm{\left(}{x}\mathit{,i}\mathrm{\right)}}\end{array}\hspace{0.17em}{\mathit{B}}_{\mathit{i}}\mathrm{\left(}{x}\mathrm{\perp }\mathit{,z}\mathrm{\right)}\mathrm{+}\mathit{n}\mathrm{\left(}{x}\mathrm{\perp }\mathrm{\right)}\mathit{.}\end{array}$(2)We assume the signal-to-noise ratio of the polarised flux to be sufficiently high for an accurate RM determination within an area A in the sky. The n-tuple of the RM values there then forms our data vector d = (RM(x _⊥ ))_{x ⊥  ∈ A}, from which all deductions will be made. The RM data can be regarded as the result of a functional linear response operator R acting on the three-dimensional magnetic field configuration and producing a two-dimensional sky image of it, $\mathit{d}\mathrm{=}\mathit{R}\hspace{0.17em}\mathit{B}\mathrm{+}\mathit{n,}$(3)where the elements of the response matrix R are defined in Eq. (2) and the application of R on B should be read as the index-sum and line-of-sight integration given in the same equation. The extended Faraday effect in the intra-cluster medium in front of the Hydra A radio galaxy was accurately measured by Taylor & Perley (1993). The intra-cluster medium has a roughly spherical symmetric electron density n_e(x _⊥ ,z) centred on the Hydra A cD galaxy. The radio jets of Hydra A inflate bubbles with relativistic, radio emitting plasma, which are visible as cavities in the X-ray emission of the surrounding hot cluster atmosphere (Wise et al. 2007). On its way to the observer the polarised radio emission of these bubbles experiences a rotation of its polarisation plane because of the Faraday effect of the magnetised intra-cluster medium.

2.2. Gaussian random fields

Owing to the information lost in the projection, Eqs. (2) and (3) are not directly invertible. It is therefore impossible to calculate the line-of-sight magnetic field B_z from φ and the known electron density profile. However, we will see that assuming statistical isotropy of the magnetic fields and using a statistical approach we are able to retrieve the amplitude spectrum of the magnetic field in Fourier space, the magnetic power spectrum $\mathit{\omega }\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\mathrm{⟨}{B}\mathrm{\left(}{k}\mathrm{\right)}\mathrm{·}\overline{{B}}\mathrm{\left(}{k}\mathrm{\right)}\mathrm{⟩}$, and the typical field strength B in position space and its autocorrelation length λ_B, which are determined by ω(k). The expectation brackets should be understood as an average over the ensemble of possible magnetic field configurations produced by typical intra-cluster dynamics. We will assume that the magnetic field statistics generated by those is Gaussian.

Gaussian magnetic fields statistics is not what numerical simulations of MHD turbulence find, see Waelkens et al. (2009). However, they are the starting point of any analysis of correlation functions. If all we know statistically about the fields is their two-point correlation, the assumption which only expresses this knowledge is a Gaussian of which this correlation matrix is the covariance matrix. Any other distribution function would contain more information in a Shannon-Boltzmann sense.

Before we can gain additional information on higher order statistics, we need to understand more the two point correlation function, and this is possible even with an approximate assumption on the statistics. Non-Gaussianities could modify the result of a power spectrum measurement based on the assumption of Gaussian fields. However, this modification can be expected to be of a modest nature, much smaller than the systematic uncertainties we deal with because of the imperfect modeling of the galaxy cluster’s atmosphere structure and the geometry of the radio source that was used to probe for Faraday rotation. Our method is relatively robust to non-Gaussian modifications, because the observational quantities it exploits, the product of RM values at different locations, are a linear transformation of the magnetic two-point correlation, irrespective of whether the underlying magnetic field statistics are Gaussian or not.

Finally, many other methods also use Gaussian random fields or assume them implicitly, e.g. Murgia et al. (2004); Bonafede et al. (2010a,b); Vacca et al. (2010); Govoni et al. (2010).

Using again the functional vector notation, we can write the Gaussian probability to find a magnetic field configuration B given a magnetic covariance matrix M =  ⟨ B B^T ⟩ , with elements M_ij(x,y) =  ⟨ B_i(x) B_j(y) ⟩ : $\mathit{P}\mathrm{\left(}\mathit{B}\mathrm{\right|}\mathit{M}\mathrm{)}\mathrm{=}\mathrm{𝒢}\mathrm{(}\mathit{B,M}\mathrm{)}\mathrm{=}\mathrm{|}\mathrm{2}\mathit{\pi }\hspace{0.17em}\mathit{M}{\mathrm{|}}^{\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}}\hspace{0.17em}\exp \left(\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}\hspace{0.17em}{\mathit{B}}^{\mathrm{T}}{\mathit{M}}^{-1}\hspace{0.17em}\mathit{B}\right)\mathit{,}$(4)with  | M |  being the determinant of the magnetic covariance matrix.

2.3. Magnetic correlation tensor

If there is a statistical translation invariance, M_ij(x,y) depends only on the distance x − y because the absolute position cannot matter and only relative positions determine the two-point statistics. As a consequence, the magnetic covariance matrix becomes diagonal in its Fourier space representation with respect to the k-vector: $\mathrm{⟨}\mathit{B̂}\mathit{i}\mathrm{\left(}{k}\mathrm{\right)}\hspace{0.17em}\overline{\mathit{B̂}\mathit{j}\mathrm{\left(}{{k}}^{\mathbf{\prime }}\mathrm{\right)}}\mathrm{⟩}\mathrm{=}\mathrm{(}\mathrm{2}\mathit{\pi }{\mathrm{)}}^{\mathrm{3}}\hspace{0.17em}\mathit{\delta }\mathrm{(}\mathit{k}\mathrm{-}{\mathit{k}}^{\mathrm{\prime }}\mathrm{)}\hspace{0.17em}\mathit{M̂}\mathit{ij}\mathrm{(}{k}\mathrm{)}\mathit{.}$(5)If we deal with statistical isotropy, one moreover finds for solenoidal fields (∇·B = 0) $\mathit{M̂}\mathit{ij}\mathrm{\left(}{k}\mathrm{\right)}\mathrm{=}\frac{\mathit{\omega }\mathrm{\left(}\mathit{k}\mathrm{\right)}}{\mathrm{2}}\hspace{0.17em}\left(\mathrm{1}\mathrm{-}\frac{{\mathit{k}}_{\mathit{i}}\hspace{0.17em}{\mathit{k}}_{\mathit{j}}}{{\mathit{k}}^{\mathrm{2}}}\right)\mathrm{+}\mathrm{i}\hspace{0.17em}{\mathit{\epsilon }}_{\mathit{ijl}}\hspace{0.17em}\frac{{\mathit{k}}_{\mathit{l}}}{\mathit{k}}\hspace{0.17em}\mathit{H}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathit{,}$(6)(see e.g. Subramanian 1999; Enßlin & Vogt 2003) and therefore $\mathit{\omega }\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\sum _{\mathit{i}}\hspace{0.17em}\mathit{M̂}\mathit{ii}\mathrm{\left(}{k}\mathrm{\right)}\mathit{.}$(7)H(k) specifies the magnetic helicity, which does not play any role in this work, because it does not lead to any detectable signature in Faraday rotation maps alone, as we will see shortly.

2.4. Faraday statistics

Assuming Gaussian statistics for turbulent magnetic fields has the pleasant property that the probability function of the Faraday depth φ is also Gaussian: $\begin{array}{ccc}\mathit{P}\mathrm{\left(}\mathit{\phi }\mathrm{\right)}& \mathrm{=}& \int \mathrm{𝒟}\mathit{B}\hspace{0.17em}\mathit{\delta }\mathrm{(}\mathit{\phi }\mathrm{-}\mathit{R}\hspace{0.17em}\mathit{B}\mathrm{)}\hspace{0.17em}\mathit{P}\mathrm{\left(}\mathit{B}\mathrm{\right)}\\ & \mathrm{=}& \int \frac{\mathrm{𝒟}\mathit{B}}{\mathrm{|}\mathrm{2}\mathit{\pi }\hspace{0.17em}\mathit{M}{\mathrm{|}}^{\frac{\mathrm{1}}{\mathrm{2}}}}\hspace{0.17em}\mathrm{\int }\mathrm{𝒟}\mathit{\eta }\exp \left(\mathrm{2}\mathit{\pi }\hspace{0.17em}\mathrm{i}\hspace{0.17em}{\mathit{\eta }}^{\mathrm{T}}\mathrm{(}\mathit{\phi }\mathrm{-}\mathit{R}\hspace{0.17em}\mathit{B}\mathrm{)}\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}\hspace{0.17em}{\mathit{B}}^{\mathrm{T}}{\mathit{M}}^{-1}\hspace{0.17em}\mathit{B}\right)\\ & \mathrm{=}& \int \mathrm{𝒟}\mathit{\eta }\exp \left(\mathrm{2}\mathit{\pi }\hspace{0.17em}\mathrm{i}\hspace{0.17em}{\mathit{\eta }}^{\mathrm{T}}\mathit{\phi }\mathrm{-}\frac{\mathrm{(}\mathrm{2}\mathit{\pi }{\mathrm{)}}^{\mathrm{2}}}{\mathrm{2}}\hspace{0.17em}{\mathit{\eta }}^{\mathrm{T}}{\mathit{R}}^{\mathrm{T}}\mathit{M}\hspace{0.17em}\mathit{R}\hspace{0.17em}\mathit{\eta }\right)\\ & \mathrm{=}& \mathrm{|}\mathrm{2}\mathit{\pi }\hspace{0.17em}{\mathit{R}}^{\mathrm{T}}\mathit{M}\hspace{0.17em}\mathit{R}{\mathrm{|}}^{\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}}\hspace{0.17em}\exp \left(\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}\hspace{0.17em}{\mathit{\phi }}^{\mathrm{T}}\mathrm{(}{\mathit{R}}^{\mathrm{T}}\mathit{M}\hspace{0.17em}\mathit{R}{\mathrm{)}}^{-1}\hspace{0.17em}\mathit{\phi }\right)\\ & \mathrm{=}& \mathrm{𝒢}\mathrm{\left(}\mathit{\phi ,}{\mathit{R}}^{\mathrm{T}}\mathit{M}\hspace{0.17em}\mathit{R}\mathrm{\right)}\mathit{.}\end{array}$(8)Here, we used the Fourier representation of the Dirac delta function (in functional spaces) as well as the well known property of Gaussian integrals , which also holds for functional or path integrals denoted by . The φ-covariance C_φ =  ⟨ φ φ^T ⟩  = R^TM R reads for our specific response and magnetic covariance matrices $\hspace{0.17em}\begin{array}{}\begin{array}{}\end{array}{\mathit{C}}_{\mathit{\phi }}\mathrm{\left(}{x}\mathrm{\perp }\mathit{,}{x}\begin{array}{}\mathrm{\prime }\\ \mathrm{\perp }\end{array}\mathrm{\right)}\mathrm{=}{\mathit{a}}_{\mathrm{0}}^{\mathrm{2}}{\mathrm{\int }}_{{\mathit{z}}_{\mathrm{1}}}^{\mathrm{\infty }}\mathrm{d}\mathit{z}{\mathrm{\int }}_{{\mathit{z}}_{\mathrm{2}}}^{\mathrm{\infty }}\mathrm{d}{\mathit{z}}^{\mathrm{\prime }}\hspace{0.17em}{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\mathrm{\right)}\hspace{0.17em}{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\mathrm{\prime }\mathrm{\right)}\hspace{0.17em}\begin{array}{}\underset{\mathrm{􏽼}0\mathrm{􏽻􏽺}0\mathrm{􏽽}}{\mathrm{⟨}\hspace{0.17em}{\mathit{B}}_{\mathit{z}}\mathrm{\left(}{x}\mathrm{\right)}{\mathit{B}}_{\mathit{z}}\mathrm{\left(}{x}\mathrm{\prime }\mathrm{\right)}\hspace{0.17em}\mathrm{⟩}}\\ {\mathit{M}}_{\mathit{zz}}\mathrm{\left(}{x}\mathit{,}{x}\mathrm{\prime }\mathrm{\right)}\end{array}\hspace{0.17em}\mathit{,}\end{array}$(9)where z₁ and z₂ denote the depth to the surface of the radio source at sky position x _⊥  and ${x}\begin{array}{}\mathrm{\prime }\\ \mathrm{\perp }\end{array}$, respectively. If the spatial electron distribution is constant, the magnetic fields are statistically homogeneous and isotropic and the observation samples a virtually infinitely extended Faraday active volume, this can be simplified even more, and reads in Fourier space $\begin{array}{ccc}\mathrm{⟨}\mathit{\phi ̂}\mathrm{\left(}{k}\mathrm{\right)}\hspace{0.17em}\overline{\mathit{\phi ̂}\mathrm{\left(}{k}\mathrm{\prime }\mathrm{\right)}}\mathrm{⟩}& \mathrm{=}& \mathrm{\left(}\mathrm{2}\mathit{\pi }{\mathrm{)}}^{\mathrm{3}}\hspace{0.17em}\mathit{\delta }\mathrm{\right(}{k}\mathrm{-}{k}\mathrm{\prime }\mathrm{\left)}\hspace{0.17em}\mathit{Ĉ}\mathit{\phi }\mathrm{\right(}\mathit{k}\mathrm{)}\\ \mathit{Ĉ}\mathit{\phi }\mathrm{\left(}\mathit{k}\mathrm{\right)}& \mathrm{=}& \frac{\mathrm{1}}{\mathrm{2}}\hspace{0.17em}{\mathit{a}}_{\mathrm{0}}^{\mathrm{2}}\hspace{0.17em}{\mathit{L}}_{\mathit{z}}\hspace{0.17em}{\mathit{n}}_{\mathrm{e}}^{\mathrm{2}}\hspace{0.17em}\mathit{\omega }\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathit{,}\end{array}$where L_z is the depth of the Faraday screen. Thus, clearly the magnetic energy, but not the helicity spectrum can be retrieved from Faraday rotation maps. In the more realistic case the geometry is less simple (finite observed volume, structured electron density, and magnetic field strength profiles, etc.). This formula has to be replaced by one which takes these effects into account, as derived in Enßlin & Vogt (2003). Because then the different Fourier modes couple, it can become simpler to again work in the position space, as we will do below (see Sect. 3.4).

2.5. REALMAF and other power spectra estimators

We model the magnetic power spectrum as a combination of spectral basis functions $\mathit{\omega }\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\sum _{\mathit{i}}{\mathit{s}}_{\mathit{i}}{\mathit{w}}_{\mathit{i}}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathit{,}$(12)where the spectral amplitudes form our spectral signal s = (s_i). The number of bins, and therefore the number of free parameters of our inference problem, will be set by a trade-off between the aim to have high spectral resolution and to have small and uncorrelated vertical error bars in the spectrum¹.

A probability function for s given the data d, the posterior P(s | d), will be set up later in this work in Sect. 3.2. Maximising the posterior with respect to the parameters in s provides an estimator for ω(k) expressed as amplitudes of relatively independent frequency bands. A similar approach has already been applied to the Hydra A cluster by Vogt & Enßlin (2005). Expressing the spectral basis functions as analogous basis functions for the magnetic autocorrelation allows us to calculate the RM covariance matrix completely in position space, which makes two rough simplifications in the processing used by Vogt & Enßlin (2005) unnecessary (see Sect. 5 for further discussion). Moreover we include noise effects into the analysis and model the Faraday screen more accurately using recently improved measurements of the geometry of the cavities in the Hydra A cluster by Wise et al. (2007).

The Hydra A cluster has been investigated by Laing et al. (2008) as well. They compared simulated rotation measure structure functions with the observed ones. A comparison with their results will be discussed in Sect. 5.

Another approach for retrieving magnetic power spectra is the one proposed by Murgia et al. (2004). Assuming a power law power spectrum, they simulated Faraday rotation maps by using the FARADAY code and compared the observed and synthetic profiles of the dispersion of the RM as fit criterion. In order to find a good fit, Murgia et al. (2004) assumed a power law to restrict the number of degrees of freedom. This method was also applied e.g. by Govoni et al. (2006) and Guidetti et al. (2008).

The REALMAF code presented in this work should alleviate several of the restrictions and assumptions of former approaches. Although REALMAF also requires the spectral range to be set up by hand, the posterior statistic calculated by it guides the choice. The power spectrum within the chosen range is relatively free to adopt any slope required by the data and permitted by the number of assumed spectral bands.

3. Inference method

3.1. Statistical inference and frequentist estimates

In this work concepts from statistical inference will be applied to Faraday rotation data to measure magnetic field properties. Because this is traditionally done via frequentist estimates, we here briefly recall the conceptual differences in the methods.

In a traditional frequentist data analysis the question is asked how likely the data would have emerged under certain model parameters, in the following called the signal. The parameters or signal values with high likelihood are then favoured over the ones with low likelihood. The likelihood is often estimated by comparing a statistical summary of the real data with that of simulated mock data sets and by characterising their similarity and difference. In our case the Faraday rotation dispersion and empirical correlation functions are typically used as statistical summaries of the RM data. A comparison of data and model seems to be transparent, because it happens completely in data space. One can easily inspect the mock data by eye and compare them to the real one. This is one reason why frequentist methods are popular. Another one is their lower computational complexity, because they usually only require calculating the measurement process in a forward fashion.

In statistical inference, a different question is addressed: How probable is a signal or a model given the data and other information sources? This is the reverse from the frequentist question and is actually the more relevant question scientifically. Answering it requires some knowledge or assumption on how probable the model was before taking the data. This prior probability, in combination with the likelihood probability, permits us – using the Bayes theorem – to construct the posterior probability of a model given the data and any earlier knowledge, which then answers the question. Because this deduction takes the reverse causal way of the measurement process, the computational complexity is usually higher than that of a comparable frequentist method.

Thus, the statistical inference (or the Bayesian) approach is more abstract in trying to do deductions in the model or signal space. Its answer is a probability function in model space and not a single estimate. The mean and variance of the posterior are often used to characterise the solution. And no mock data and statistical summaries are usually generated, because the inference inserts the observed data into a likelihood function from which a backward deduction of the model is done.

Frequentist and Bayesian results often agree. However, there are at least three reasons why they can differ, and we believe that the latter two matter here:

• 1.

  Bayesian methods specify the prior assumptions explicitely,while those are often implicit in frequentist estimates. Becausethe prior has some impact on the result, different choices of priorswill lead to different results.

• 2.

  Bayesian methods are better suited to exploit the full information content of the data. The summary statistics used in frequentist analysis often looses a significant amount of the information imprinted in the data. In our case, we are investigating the correlation between all pairs of pixels, but a typical frequentist statistics might only look at the variance of the individual pixels or the correlation function from a subset of similar pixels, and therefore exploits much less of the data complexity.

• 3.

  There are situations where no suitable summary statistics exist to construct an empirical likelihood. In our case, calculating an empirical RM correlation function to compare data with mock data only makes sense if the RM field can be assumed to result from a statistically homogeneous process. However, we are facing a very inhomogeneous geometry here because of the highly structured cluster atmosphere, with very few comparable lines of sight. How to obtain a meaningful two-point correlation function from the data in this situation is not clear to us.

Because of the last point, the Bayesian inference performed in this work does not provide common quantities of frequentist analysis like mock data, observed RM correlation function, and the like. The RM data are confronted with the theoretical covariance matrix that is expected for a given power spectrum of the magnetic field after a model of the inhomogeneous cluster atmosphere is factored out. The power spectrum would, in a statistically homogeneous setting, describe an easily observable RM correlation function, but not in our inhomogeneous case.

The Bayesian approach also extracts statistical information from poorly sampled quantities, while providing the correct large uncertainty resulting from the undersampling. For example in our case the magnetic modes on scales of or larger than the observed RM map are clearly undersampled. The largest modes appear as a single number in the data, the mean RM value of the map. Because it is a single number, it is hard to invert this into precise spectral information on large-scale magnetic fluctuations. Nevertheless, there is some useful information in this number, especially if it is far from zero. In this case, the total magnetic spectral power on large spatial scales should be on a level which makes the observed mean RM value sufficiently likely. The use of a Jeffreys prior on all spectral amplitudes enforces in our application that the inversion of poorly constrained modes is made in a very conservative fashion. With frequentist methods the usage of the information in such poorly sampled modes is more difficult, because no empirical statistics can be built up from a single observation.

3.2. Statistical basics

Our aim is to retrieve from the rotation measures (data d) the magnetic power spectrum represented by the spectral amplitudes (signal s). The Bayes theorem expresses the posterior P(s | d) (the conditional probability of the signal given the data) in terms of the likelihood P(d | s) in the following way: $\mathit{P}\mathrm{\left(}\mathit{s}\mathrm{\right|}\mathit{d}\mathrm{)}\mathrm{=}\frac{\mathit{P}\mathrm{\left(}\mathit{d}\mathrm{\right|}\mathit{s}\mathrm{)}\hspace{0.17em}\mathit{P}\mathrm{(}\mathit{s}\mathrm{)}}{\mathit{P}\mathrm{\left(}\mathit{d}\mathrm{\right)}}\mathrm{·}$(13)Because P(d) is a constant with respect to s, the most probable s given d can be retrieved by maximising the joint probability, which consists of the likelihood P(d | s) multiplied with the prior P(s). These will be specified now for our problem.

In our case the signal is the ensemble of spectral amplitudes s = (s_i)_i. We adopt Jeffreys prior for each spectral amplitude, so that $\mathit{P}\mathrm{\left(}\mathit{s}\mathrm{\right)}\mathrm{\propto }{{}^{\mathrm{􏽑}}}_{\mathit{i}}{\mathit{s}}_{\mathit{i}}^{-1}$. Jeffreys prior is a uniform distribution on a logarithmic scale, expressing that we are beforehand not even sure about the magnitude of the magnetic spectrum.

Our likelihood P(d | s) has to incorporate two stochastic processes. The first is the realisation of the magnetic field configuration according to Eq. (4). Because s fully determines the here relevant symmetric part of the magnetic correlation tensor, M = M(ω,H) = M(∑ _is_i ω_i,H), the statistics of the Faraday depth is determined by it as well, and it is also Gaussian according to Eq. (8). If the noise is also Gaussian, then the data are Gaussian with a covariance matrix C, with entries C_ij =  ⟨ d_i d_j ⟩  for the correlation of RM data at locations x _⊥ ,i and x _⊥ ,j. Using vector and matrix notation for the data vector d = φ + n and its covariance matrix C, we can write the Gaussian as $\mathit{P}\mathrm{\left(}\mathit{d}\mathrm{\right|}\mathit{s}\mathrm{)}\mathrm{=}\mathrm{𝒢}\mathrm{(}\mathit{d,C}\mathrm{\left(}\mathit{s}\mathrm{\right)}\mathrm{)}\mathrm{=}\mathrm{|}\mathrm{2}\mathit{\pi }\hspace{0.17em}\mathit{C}{\mathrm{|}}^{\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}}\hspace{0.17em}{\mathrm{e}}^{\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}{\mathit{d}}^{\mathit{T}}{\mathit{C}}^{-1}\mathit{d}}\mathit{.}$(14)An expression for C as a function of our spectral parameter s is needed. In the view of the linearity of Eq. (2), the observed data RM can be expressed in terms of the magnetic field B (vector over physical space), a response matrix R (projecting from physical to data space) and noise n (in data space) as: d = R B + n. From this we find $\begin{array}{ccc}\mathit{C}\mathrm{\left(}\mathit{s}\mathrm{\right)}& \mathrm{=}& \mathrm{⟨}\mathrm{(}\mathit{R}\hspace{0.17em}\mathit{B}\mathrm{+}\mathit{n}\mathrm{)}\hspace{0.17em}\mathrm{(}\mathit{R}\hspace{0.17em}\mathit{B}\mathrm{+}\mathit{n}{\mathrm{)}}^{\mathit{T}}\mathrm{⟩}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\mathit{R}\hspace{0.17em}\mathrm{⟨}\mathit{B}\hspace{0.17em}{\mathit{B}}^{\mathrm{T}}\mathrm{⟩}\hspace{0.17em}{\mathit{R}}^{\mathrm{T}}\mathrm{+}\mathrm{⟨}\mathit{n}{\mathit{n}}^{\mathrm{T}}\mathrm{⟩}\\ & \mathrm{=}& \mathit{R}\hspace{0.17em}\mathit{M}\hspace{0.17em}{\mathit{R}}^{†}\mathrm{+}\mathit{N}\mathrm{=}\hspace{0.17em}{\mathit{C}}_{\mathit{\phi }}\mathrm{\left(}\mathit{s}\mathrm{\right)}\mathrm{+}\mathit{N}\mathit{.}\end{array}$Thus to any derived Faraday rotation covariance matrix C_φ =  ⟨ φ φ^T ⟩  an estimate for the noise covariance matrix N =  ⟨ nn^T ⟩  needs to be added.

3.3. Maximising the posterior

The power spectrum inference within the maximum-a-posteriori approach can only be conducted numerically. Here, we maximise the logarithm of the posterior, because there are convenient analytical expressions for the gradient and Hessian. Multiplication of the likelihood with the prior to get the posterior is an addition in logarithmic space, which also simplifies matters. Indices i,j refer to the degrees of freedom of the parameterisation of s. We find $\frac{\mathit{\partial }\ln \mathrm{\left[}\mathit{P}\mathrm{\right(}\mathit{s}\mathrm{\left|}\mathit{d}\mathrm{\right)}\mathrm{]}}{\mathit{\partial }{\mathit{s}}_{\mathit{i}}}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}\mathrm{Tr}\left[\mathrm{(}\mathit{d}{\mathit{d}}^{\mathrm{T}}\mathrm{-}{\mathit{C}}^{\mathrm{)}}\hspace{0.17em}\left({\mathit{C}}^{-1}\frac{\mathit{\partial C}}{\mathit{\partial }{\mathit{s}}_{\mathit{i}}}{\mathit{C}}^{-1}\right)\right]\mathrm{-}\frac{\mathrm{1}}{{\mathit{s}}_{\mathit{i}}}\mathit{,}$(15)and if s contributes only linearly to C: $\begin{array}{ccc}\frac{{\mathit{\partial }}^{\mathrm{2}}\ln \mathrm{\left[}\mathit{P}\mathrm{\right(}\mathit{s}\mathrm{\left|}\mathit{d}\mathrm{\right)}\mathrm{]}}{\mathit{\partial }{\mathit{s}}_{\mathit{j}}\mathit{\partial }{\mathit{s}}_{\mathit{i}}}& \mathrm{=}& \mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}\mathrm{Tr}\left[{\mathit{C}}^{-1}\frac{\mathit{\partial C}}{\mathit{\partial }{\mathit{s}}_{\mathit{j}}}{\mathit{C}}^{-1}\frac{\mathit{\partial C}}{\mathit{\partial }{\mathit{s}}_{\mathit{i}}}\right]\\ & & \mathrm{-}\mathrm{Tr}\left[\mathrm{(}\mathit{d}{\mathit{d}}^{\mathrm{T}}\mathrm{-}{\mathit{C}}^{\mathrm{)}}\left({\mathit{C}}^{-1}\frac{\mathit{\partial C}}{\mathit{\partial }{\mathit{s}}_{\mathit{i}}}{\mathit{C}}^{-1}\frac{\mathit{\partial C}}{\mathit{\partial }{\mathit{s}}_{\mathit{j}}}{\mathit{C}}^{-1}\right)\right]\mathrm{+}\frac{{\mathit{\delta }}_{\mathit{ij}}}{{\mathit{s}}_{\mathit{i}}^{\mathrm{2}}}\mathrm{·}\end{array}$(16)The most probable power spectrum is then found iteratively. We use a second-order Newton method to solve for it, with

${\mathit{s}}_{\mathrm{new}}\mathrm{=}\mathit{s}\mathrm{-}{\left(\frac{{\mathit{\partial }}^{\mathrm{2}}\ln \mathrm{\left[}\mathit{P}\mathrm{\right(}\mathit{s}\mathrm{\left|}\mathit{d}\mathrm{\right)}\mathrm{]}}{\mathit{\partial }{\mathit{s}}_{\mathit{j}}\mathit{\partial }{\mathit{s}}_{\mathit{i}}}\right)}^{-1}\frac{\mathit{\partial }\ln \mathrm{\left[}\mathit{P}\mathrm{\right(}\mathit{s}\mathrm{\left|}\mathit{d}\mathrm{\right)}\mathrm{]}}{\mathit{\partial }{\mathit{s}}_{\mathit{i}}}\mathit{,}$where repeated indices are assumed to be summed over.

3.4. The covariance matrix

An expression for the covariance matrix ${\mathit{C}}_{\mathit{\phi }}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\mathrm{⟨}\mathit{\phi }\mathrm{\left(}{x}{}_{\mathrm{\perp }}\mathrm{\right)}\hspace{0.17em}\mathit{\phi }\mathrm{\left(}{x}\begin{array}{}\mathrm{\prime }\\ \mathrm{\perp }\end{array}\mathrm{\right)}\mathrm{⟩}$ of the rotation measure map needs to be derived in the following form: ${\mathit{C}}_{\mathit{\phi }}\mathrm{(}{x}\mathrm{\perp }\mathit{,}{x}\begin{array}{}\mathrm{\prime }\\ \mathrm{\perp }\end{array}\mathrm{)}\mathrm{=}{\mathit{a}}_{\mathrm{0}}^{\mathrm{2}}\hspace{0.17em}⟨\hspace{0.17em}{\mathrm{\int }}_{{\mathit{z}}_{\mathrm{1}}}^{\mathrm{\infty }}\mathrm{d}\mathit{z}\hspace{0.17em}{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\mathrm{\right)}{\mathit{B}}_{\mathit{z}}\mathrm{\left(}{x}\mathrm{\right)}{\mathrm{\int }}_{{\mathit{z}}_{\mathrm{2}}}^{\mathrm{\infty }}\mathrm{d}{\mathit{z}}^{\mathrm{\prime }}\hspace{0.17em}{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\mathrm{\prime }\mathrm{\right)}{\mathit{B}}_{\mathit{z}}\mathrm{\left(}{x}\mathrm{\prime }\mathrm{\right)}\hspace{0.17em}⟩\mathit{.}$The electron density distribution n_e(x) needs to be obtained from independent data, e.g. from X-ray observations. The magnetic field energy is assumed to follow n_e(x) according to Dolag et al. (2002): $\mathrm{⟨}{\mathit{B}}^{\mathrm{2}}\mathrm{\left(}{x}\mathrm{\right)}\mathrm{⟩}\mathrm{=}{\mathit{B}}_{\mathrm{0}}^{\mathrm{2}}\hspace{0.17em}\mathrm{\left(}{\mathit{n}}_{\mathrm{e}}\mathrm{\right(}{x}\mathrm{)}\mathit{/}{\mathit{n}}_{\mathrm{0}}{\mathrm{)}}^{\mathrm{2}\mathit{\alpha }}\mathit{.}$(17)The following approach is based on statistical homogeneity of the magnetic field. Therefore we re-scale  ⟨ B²(x) ⟩  with a rescaleing function h, h(x) =  ⟨ B²(x) ⟩ ^0.5/B₀ = (n_e(x)/n₀)^α. This ensures $\mathrm{⟨}\begin{array}{}{}_{\mathrm{􏽥}}\\ {\mathit{B}}^{\mathrm{2}}\end{array}\mathrm{\left(}\mathit{x}\mathrm{\right)}\mathrm{⟩}\mathrm{=}\mathrm{⟨}{\mathit{B}}^{\mathrm{2}}\mathrm{\left(}{x}\mathrm{\right)}\mathrm{⟩}\mathit{/}{\mathit{h}}^{\mathrm{2}}\mathrm{\left(}{x}\mathrm{\right)}\mathrm{=}{\mathit{B}}_{\mathrm{0}}^{\mathrm{2}}$ to be independent of position. If we also assume statistical isotropy, we obtain the scaling of the required z-component of B: $\mathrm{⟨}{\mathit{B}}_{\mathit{z}}^{\mathrm{2}}\mathrm{\left(}{x}\mathrm{\right)}\mathrm{⟩}\mathit{/}{\mathit{B}}_{\mathit{z}\mathrm{0}}^{\mathrm{2}}\mathrm{=}\mathrm{⟨}{\mathit{B}}^{\mathrm{2}}\mathrm{\left(}{x}\mathrm{\right)}\mathrm{⟩}\mathit{/}{\mathit{B}}_{\mathrm{0}}^{\mathrm{2}}$. Thus, $\begin{array}{ccc}{\mathit{C}}_{\mathit{\phi }}& \mathrm{=}& {\mathit{a}}_{\mathrm{0}}^{\mathrm{2}}{\mathrm{\int }}_{{\mathit{z}}_{\mathrm{1}}}^{\mathrm{\infty }}\mathrm{d}\mathit{z}\hspace{0.17em}{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\mathrm{\perp }\mathit{,z}\mathrm{\right)}\hspace{0.17em}\mathit{h}\mathrm{\left(}{x}\mathrm{\perp }\mathit{,z}\mathrm{\right)}\\ & & \int \begin{array}{}\mathrm{\infty }\\ {\mathit{z}}_{\mathrm{2}}\end{array}\mathrm{d}{\mathit{z}}^{\mathrm{\prime }}\hspace{0.17em}{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\begin{array}{}\mathrm{\prime }\\ \mathrm{\perp }\end{array}\mathit{,}{\mathit{z}}^{\mathrm{\prime }}\mathrm{\right)}\mathit{h}\mathrm{\left(}{x}\begin{array}{}\mathrm{\prime }\\ \mathrm{\perp }\end{array}\mathit{,}{\mathit{z}}^{\mathrm{\prime }}\mathrm{\right)}\mathrm{⟨}\begin{array}{}\mathrm{􏽥}\\ {\mathit{B}}_{\mathit{z}}\end{array}\mathrm{\left(}{x}\mathrm{\perp }\mathit{,z}\mathrm{\right)}\begin{array}{}\mathrm{􏽥}\\ {\mathit{B}}_{\mathit{z}}\end{array}\mathrm{\left(}{x}\begin{array}{}\mathrm{\prime }\\ \mathrm{\perp }\end{array}\mathit{,}{\mathit{z}}^{\mathrm{\prime }}\mathrm{\right)}\mathrm{⟩}\\ & \mathrm{=}& {\mathit{a}}_{\mathrm{0}}^{\mathrm{2}}{\mathit{n}}_{\mathrm{0}}^{\mathrm{2}}{\mathrm{\int }}_{{\mathit{z}}_{\mathrm{1}}}^{\mathrm{\infty }}\mathrm{d}\mathit{z}\hspace{0.17em}\mathit{f}\mathrm{\left(}{x}\mathrm{\perp }\mathit{,z}\mathrm{\right)}{\mathrm{\int }}_{{\mathit{z}}_{\mathrm{2}}\mathrm{-}\mathit{z}}^{\mathrm{\infty }}\mathrm{d}{\mathit{r}}_{\mathit{z}}\hspace{0.17em}\mathit{f}\mathrm{(}{x}\begin{array}{}\mathrm{\prime }\\ \mathrm{\perp }\end{array}\mathit{,}{\mathit{r}}_{\mathit{z}}\mathrm{+}\mathit{z}\mathrm{)}\hspace{0.17em}{\mathit{M}}_{\mathit{zz}}\mathrm{\left(}{r}\mathrm{\right)}\mathit{,}\end{array}$(18)with ${r}\mathrm{=}\mathrm{(}{x}{}_{\mathrm{\perp }}\mathrm{-}{x}\begin{array}{}\mathrm{\prime }\\ \mathrm{\perp }\end{array}\mathit{,}{\mathit{r}}_{\mathit{z}}\mathrm{)}$, where we introduced the window function f(x), which is nonzero inside the Faraday screen: $\mathit{f}\mathrm{\left(}{x}\mathrm{\right)}\mathrm{=}\frac{{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\mathrm{\right)}}{{\mathit{n}}_{\mathrm{0}}}\mathit{h}\mathrm{\left(}{x}\mathrm{\right)}\mathrm{=}{\left(\frac{{\mathit{n}}_{\mathrm{e}}\mathrm{\left(}{x}\mathrm{\right)}}{{\mathit{n}}_{\mathrm{0}}}\right)}^{\mathit{\alpha }\mathrm{+}\mathrm{1}}\mathrm{·}$(19)${\mathit{M}}_{\mathit{zz}}\mathrm{\left(}{r}\mathrm{\right)}\mathrm{=}\mathrm{⟨}\begin{array}{}{}_{\mathrm{􏽥}}\\ {\mathit{B}}_{\mathit{z}}\end{array}\mathrm{\left(}{x}\mathrm{\right)}\begin{array}{}{}_{\mathrm{􏽥}}\\ {\mathit{B}}_{\mathit{z}}\end{array}\mathrm{(}{x}\mathrm{+}{r}\mathrm{)}\mathrm{⟩}$ is the zz-part of the (rescaled) magnetic autocorrelation tensor, which we try to measure.

3.5. Modelling the magnetic autocorrelation

M_zz(r) can be split up into a longitudinal part M_L(r) and a corresponding perpendicular part M_N(r), which depend only on the absolute value of r (see Subramanian 1999; or Enßlin & Vogt 2003). ${\mathit{M}}_{\mathit{zz}}\mathrm{\left(}{r}\mathrm{\right)}\mathrm{=}{\mathit{M}}_{\mathrm{L}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\frac{{\mathit{r}}_{\mathit{z}}^{\mathrm{2}}}{{\mathit{r}}^{\mathrm{2}}}\mathrm{+}{\mathit{M}}_{\mathrm{N}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\frac{{\mathit{r}}_{\mathrm{\perp }}^{\mathrm{2}}}{{\mathit{r}}^{\mathrm{2}}}\mathrm{·}$(20)The rescaled magnetic autocorrelation function $\mathit{\omega }\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\mathrm{⟨}\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{\left(}{x}\mathrm{\right)}\mathrm{·}\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{(}{x}\mathrm{+}{r}\mathrm{)}\mathrm{⟩}$ is given by $\mathit{\omega }\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\sum _{\mathit{i}}{\mathit{M}}_{\mathit{ii}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\mathrm{2}{\mathit{M}}_{\mathrm{N}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{+}{\mathit{M}}_{\mathrm{L}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathit{.}$(21)The information of the data can be combined into a one-dimensional correlation function and therefore does not allow us to simultaneously determine two independent functions. Thus a relation between M_N and M_L is necessary. There are two suitable possibilities how to connect them.

• a)

  Full isotropy assumption, which additionally assumesisotropy at any location in Fourier space. This yields${\mathit{M}}_{\mathrm{N}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}{\mathit{M}}_{\mathrm{L}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathit{.}$(22)

• b)

  Considering the solenoidal character of the magnetic field by $\mathrm{\nabla }\mathrm{·}\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{=}\mathrm{0}$, which yields ${\mathit{M}}_{\mathrm{N}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}\mathit{r}}\frac{\mathrm{d}}{\mathrm{d}\mathit{r}}\mathrm{\left(}{\mathit{r}}^{\mathrm{2}}{\mathit{M}}_{\mathrm{L}}\mathrm{\right(}\mathit{r}\mathrm{\left)}\mathrm{\right)}\mathit{.}$(23)

The $\mathrm{\nabla }\mathrm{·}\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{=}\mathrm{0}$ condition leads to Eq. (6) and is more physical, but also more complicated to handle. Both assumptions have been implemented in REALMAF. However, the most computationally efficient way is to use the full isotropy assumption and apply the divergence relation in Fourier space, which results in a simple multiplication factor, see Sect. 3.6.

Connecting of M_L and M_N allows us to parameterise M_zz(r) using ${\mathit{M}}_{\mathit{zz}}\mathrm{\left(}{r}\mathrm{\right)}\mathrm{=}\sum _{\mathit{i}}{\mathit{s}}_{\mathit{i}}{\mathit{M}}_{\mathit{zz}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}{r}\mathrm{\right)}\mathit{,}$(24)where ${\mathit{M}}_{\mathit{zz}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}{r}\mathrm{\right)}$ is further split into ${\mathit{M}}_{\mathrm{L}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{r}\mathrm{\right)}$ and ${\mathit{M}}_{\mathrm{N}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{r}\mathrm{\right)}$ described by Eq. (20). Because of the linearity of the Fourier transformation, the linear coefficients s_i in real space also appear as linear coefficients in the analogous representation in Fourier space. Thus, although we are calculating the covariance matrix completely in real space, the power spectrum can be obtained directly:

$\mathit{\omega }\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\sum _{\mathit{i}}\hspace{0.17em}{\mathit{s}}_{\mathit{i}}\hspace{0.17em}\mathrm{\left(}\mathrm{2}{\mathit{M}}_{\mathrm{N}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\right(}\mathit{k}\mathrm{)}\mathrm{+}{\mathit{M}}_{\mathrm{L}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{(}\mathit{k}\mathrm{\left)}\mathrm{\right)}\mathrm{=}\sum _{\mathit{i}}\hspace{0.17em}{\mathit{s}}_{\mathit{i}}\hspace{0.17em}{\mathit{w}}_{\mathit{i}}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathit{.}$The basis functions together with the retrieved spectral amplitudes s_i can be translated into ω(r) (see Fig. 1a) and the corresponding ω(k) (see Fig. 1b). The used model functions for ${\mathit{M}}_{\mathrm{N}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}$ and ${\mathit{M}}_{\mathrm{L}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}$ were chosen to transform into relatively well separated and well localised frequency bands w_i(k) in Fourier space. Appendix A shows the used model functions and describes how they were derived.

3.6. The divergence relation in Fourier space

The zz-part of the magnetic autocorrelation tensor can be expressed in Fourier space analogously to Eq. (20) ${\mathit{M}}_{\mathit{zz}}\mathrm{\left(}{k}\mathrm{\right)}\mathrm{=}\begin{array}{}\mathrm{􏽢}\\ {\mathit{M}}_{\mathrm{L}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}\frac{{\mathit{k}}_{\mathit{z}}^{\mathrm{2}}}{{\mathit{k}}^{\mathrm{2}}}\mathrm{+}\begin{array}{}\mathrm{􏽢}\\ {\mathit{M}}_{\mathrm{N}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}\frac{{\mathit{k}}_{\mathrm{\perp }}^{\mathrm{2}}}{{\mathit{k}}^{\mathrm{2}}}\mathrm{·}$(25)Applying $\mathrm{\nabla }\mathrm{·}\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{=}\mathrm{0}$ gives then $\begin{array}{}\mathrm{􏽢}\\ {\mathit{M}}_{\mathrm{L}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\mathrm{0}\hspace{0.17em}\mathrm{=}\mathrm{\Rightarrow }\hspace{0.17em}\begin{array}{}\mathrm{􏽢}\\ \mathit{\omega }\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\mathrm{2}\begin{array}{}\mathrm{􏽢}\\ {\mathit{M}}_{\mathrm{N}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathit{,}$(26)see Enßlin & Vogt (2003). The $\begin{array}{}{}_{\mathrm{􏽢}}\\ {\mathit{M}}_{\mathrm{L}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}$ from Eq. (25) is in general not the Fourier transformation of M_L(r) used in Eq. (20). To create a connection to the real-space parameterisation M_zz(r), we multiply $\begin{array}{}{}^{\mathrm{􏽢}}\\ \mathit{\omega }\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}$ by 3/2, which is analogous to an addition of $\begin{array}{}{}_{\mathrm{􏽢}}\\ {\mathit{M}}_{\mathrm{L}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}$ with $\begin{array}{}{}_{\mathrm{􏽢}}\\ {\mathit{M}}_{\mathrm{L}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\begin{array}{}{}_{\mathrm{􏽢}}\\ {\mathit{M}}_{\mathrm{N}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}$. This makes M_zz(k) a spherically symmetric function, which forces that its Fourier back-transform M_zz(r) must be spherically symmetric as well. Then, ${\mathit{M}}_{\mathrm{L}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}{\mathit{M}}_{\mathrm{N}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\hspace{0.17em}\mathrm{\Leftarrow }\mathrm{\Rightarrow }\hspace{0.17em}\begin{array}{}\mathrm{􏽢}\\ {\mathit{M}}_{\mathrm{L}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\begin{array}{}\mathrm{􏽢}\\ {\mathit{M}}_{\mathrm{N}}\end{array}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathit{,}$(27)which is the full isotropy assumption defined in Eq. (22). Therefore, when assuming full isotropy, the divergence condition can be imposed on ω(k) to a good approximation via multiplying it with the factor 2/3. This was verified as well by numerical experiments². Using the full isotropy assumption reduces the processing time by up to a factor 4 and allows a simpler spherically symmetric model for M_zz. In practice the power spectrum is multiplied with a factor 2/3 and the resulting magnetic field strength by $\sqrt{\mathrm{2}\mathit{/}\mathrm{3}}$.

Finally, the $\mathrm{\nabla }\mathrm{·}\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{=}\mathrm{0}$ condition, which is used by Vogt & Enßlin (2005) as well, is not exact because of

$\mathrm{\nabla }\mathrm{·}\begin{array}{}\mathrm{􏽥}\\ \mathit{B}\end{array}\mathrm{=}\mathit{h}\mathrm{\nabla }\mathrm{·}\mathit{B}\mathrm{+}\mathit{B}\mathrm{·}\mathrm{\nabla }\mathit{h}\mathrm{=}\mathit{B}\mathrm{·}\mathrm{\nabla }\mathit{h,}$with the rescaleing function h defined in Sect. 3.4. If B varies on scales much smaller than the ones on which h varies, this assumption is sufficiently close to the condition ∇·B = 0. Another indication that the resulting errors are negligibl small is that a totally different assumption on the statistical distribution of the magnetic fields, like full isotropy, which violates $\mathrm{\nabla }\mathrm{·}\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{=}\mathrm{0}$ by definition, produces only a shift by a factor of 1.5, as shown above.

4. Application to the Hydra A data

4.1. The data

High-resolution polarisation measurements of the radio lobes of Hydra A at different wavelengths, with the VLA suitable for RM map making, were performed by Taylor & Perley (1993). High-fidelity RM maps were constructed from these using the PACERMAN code (Dolag et al. 2005; Vogt et al. 2005). PACERMAN uses error weighting and also non-local information to solve the nπ-ambiguity and therefore produces more reliable Faraday rotation maps than algorithms that do standard fitting of Eq. (2) pixel by pixel. The resulting RM map was tested by Vogt et al. (2005) for the presence of map-making artifacts. These artifacts can be statistically detected using the anti-correlation of inferred RM and position angle values imprinted by misfits, as first shown by Enßlin et al. (2003). The tests show that the current RM map of the south lobe of Hydra A is prone to nπ-ambiguities and other artifacts, even if generated with PACERMAN. However, the north lobe, which contains lower RM values and is therefore pointing towards us, is sufficiently clean. Therefore Vogt & Enßlin (2005) analysed only the north lobe of Hydra A, and our application of REALMAF is also restricted to this area.

The resolution of the RM map is too high to take every single pixel into account in our analysis. A computational bottleneck is the double integral in the calculation of the model covariance matrix followed by matrix multiplications described in Sect. 3.3. To reduce the amount of required computational resources, we merge areas in the map to bigger pixels by combining the rotation measure values, the errors, and the pixel position in an inverse noise-weighted fashion. The idea is to average the data in quadratic areas chosen so that they have roughly the same noise level. Therefore the outer regions with lower radio surface brightness and consequently larger RM errors are more heavily averaged. For details on the algorithm combining the pixels, the reader is referred to Vogt & Enßlin (2005).

These merged data points then form the data vector d used in Sect. 3.2. Figure 2 visualises the merged map by showing the averaged values in colour within the square areas from which they were averaged. The area outside the squares does not contain any used RM values. In Fig. 3 the radial profile of of RM² is shown and is compared with the ensemble averaged  ⟨ RM² ⟩  predicted by the deduced magnetic spectrum and assumed radial profile of the electron and magnetic field energy densities. From this figure it can be seen that the largest separation between our pixels is r_max = 33 kpc, whereas the smallest r_min = 0.073 kpc. A first guess of the optimal spectral range can be obtained based on these maximal and minimal sizes of visible structures calculated using k_min = π/r_max and k_max = π/r_min The final result in Fig. 7 uses roughly these spectral boundaries.

Spatial averaging of RM data is not properly modelled in Eq. (18). In principle, it could be included into this formula and the generic inference formalism, owing to the linearity of averaging. However, we have to perform the double integral in Eq. (18) numerically for any pair of data pixels. Adding four more integration dimensions would exceed our computational possibilities.

We accept the spectral error resulting from the data averaging operation because pixels that resulted from averaging over large regions are automatically only contributing information about scales above their size. Any baseline formed by them with other pixels will be of the order of the pixel size or larger. The product of the RM values of a pair of pixels is compared with the expected data correlation over this baseline distance. Thus these pixels determine the inferred spectrum mostly on scales on which their internal RM substructure, which was averaged out, does not matter.

The silhouette of the radio lobe used to probe the Faraday screen breaks the translational invariance of our problem, in addition to the effect of the electron density profile. However, this does not impact directly our data covariance matrix C, because this only contains the expected (averaged over an ensemble of data realisations) correlation of individual pairs of positions in the RM map. The silhouette restricts the available positions and therefore the number of baselines, which means that on larger scales, we have less information. However, the Bayesian formalism is informed about this and just increases the uncertainty margins for poorly constrained data modes.

We also require to specify the noise covariance matrix. Assuming the noise at different pixels to be uncorrelated, the noise matrix is diagonal. The diagonal elements can be approximated by the square of the RM-deviation of each pixel given by the error map produced by PACERMAN.

The used rotation measure map has a mean value of about 600 rad/m². The galactic rotation measure at the position of Hydra A is about 1 percent of this and therefore negligible, see Johnston-Hollitt et al. (2004). Therefore we assume that the mean value is caused by cosmic variance which is in turn processing a too small sample of the cluster. This mean value is correctly interpreted by our method in that it adds power to the small k part of the magnetic power spectrum, compared with the case in which we subtract the mean value from the data d.

4.2. Window function

The window function f defined in Eq. (19) requires the electron density distribution n_e(x). Hydra A is a cooling-core cluster. Therefore we assume a spherically symmetric double beta profile defined as ${\mathit{n}}_{\mathrm{e}}\mathrm{\left(}\mathit{x}\mathrm{\right)}\mathrm{=}{\left[{\mathit{n}}_{\mathrm{e}\mathrm{1}}^{\mathrm{2}}{\left(\mathrm{1}\mathrm{+}{\left(\frac{\mathit{r}}{{\mathit{r}}_{\mathrm{c}\mathrm{1}}}\right)}^{\mathrm{2}}\right)}^{\mathrm{-}\mathrm{3}\mathit{\beta }}\mathrm{+}{\mathit{n}}_{\mathrm{e}\mathrm{2}}^{\mathrm{2}}{\left(\mathrm{1}\mathrm{+}{\left(\frac{\mathit{r}}{{\mathit{r}}_{\mathrm{c}\mathrm{2}}}\right)}^{\mathrm{2}}\right)}^{\mathrm{-}\mathrm{3}\mathit{\beta }}\right]}^{\mathrm{1}\mathit{/}\mathrm{2}}\mathit{,}$(28)with n_e1 = 0.056 cm^-3, n_e2 = 0.0063 cm^-3, r_c1 = 33.3 kpc, r_c2 = 169 kpc and β = 0.766, as fitted by Vogt & Enßlin (2005) using data from Mohr et al. (1999). Various works propose that the observed radio lobes of active galactic nuclei, like Hydra A, form cavities in the cluster atmosphere, see e.g. Böhringer et al. (1993), Enßlin (1999), McNamara et al. (2000), Fabian et al. (2000), and David et al. (2001). In these cavities the polarised radio light is produced, but n_e(x) is presumably so low that no significant Faraday rotation is expected to occur there. The fact n_e(x) = 0 for any x ∈  cavity can be taken into account by adjusting the integration boundaries in Eq. (18). In order to calculate the important lower integration boundary, where the Faraday active area starts, we need to know the geometry of the cavity. The cavity in our used map-section can be approximated by an ellipse with the projected distance 24.9 kpc from cluster centre, projected semimajor axis 20.5 kpc, and semiminor axis with 12.4 kpc, as Wise et al. (2007) propose. Additionally the Hydra A north lobe is tilted by an angle θ towards the observer. Depending on θ the innermost radius of the probed Faraday active area shifts, as can be seen in Fig. 4. There it becomes apparent that the sensitivity of our analysis starts at radii of about 20 to 30 kpc from the cluster centre depending on θ. θ is not well determined and its estimates range from about 30° (Taylor & Perley (1993), if the cavity effect is added) to more than 45° (Laing et al. 2008) to 60° (Wise et al. 2007).

4.3. Measuring the magnetic profile

The index α of the magnetic field can be retrieved from the data. In contrast to the power spectrum³, which is built dynamically by our algorithm, α has to be set up for every processing run. Maximising the posterior with respect to α provides indications on most probable configurations. One drawback is that α also depends on the low-frequency cut-off k_min (position of the top of the first spectral bin) chosen; meaning a two-dimensional problem needs to be optimised. The result we find is α = 1 with k_min = 0.42 kpc^-1. Because this result only depends on the large scales, a rather low-resolution study with 300 data points was sufficient and permitted a fast investigation of the 2-d parameter space. Sampling of k_min for three different projection angles θ can be found in Fig. 5.

Cut-offs with k < 0.42 kpc^-1 correspond to α > 1 and vice versa. Smaller θ implies a smaller physical distance between the data points, which results in a small shift of k_min. However, it is negligible and inside the statistical variance of about 0.05 kpc^-1. Figure 6 shows a consistency check by sampling α for the given k_min = 0.42 kpc^-1. The curvature in Fig. 6 gives roughly a 1-sigma variance for α of about Δα =  ± 0.3, while assuming Gaussianity around the peak. When we also consider the impact of the variance of k_min, we find α = 1 ± 0.5.

Surprisingly the posterior maximum is rising for rising θ with no apparent limit. We suppose that this has nothing to do with the probability of the physical model and therefore θ cannot be restricted using the data of the Hydra A north lobe only. Indeed, θ is a parameter, that only affects the integration bounds of Eq. (18), whereas the power spectrum, the magnetic scaling α, and the cut-offs affect the integrand.

4.4. Power spectrum and magnetic characteristics

The shape of the power spectrum depends on the spectral cut-offs k_min and k_max. We retrieved the most probable k_min to be 0.42 kpc^-1 in Sect. 4.3. Obviously, the most probable configuration is a power law with no large-scale turnover visible within the probed spatial scales, as can be seen in Fig. 7. The decline for smaller k_min may be caused by the limited rotation measure map dimensions and the down-pull of thereby unconstrained amplitudes by Jeffreys prior. A larger k_min will remove necessary degrees of freedom of the model at the large scales, which forces the algorithm to deposit additional power at the neighbouring scale. The high-frequency cut-off is limited by the resolution of the map and instrumental beam. A beam with spatial scales of 0.32 kpc (FWHM) gives k_max = 10 kpc^-1. Below, we increase the resolution by increasing the number of data points to 1000 and extend the spectrum to a maximal possible range. We made three processing runs with α = 1 and three possible projection angles, see Fig. 8 and Table 1. To investigate the dependence of the results on α we made two more runs with α = 0.5 and α = 1.5, see Table 2. Formulas for the calculation of the magnetic characteristics are described in Appendix B. The given statistical uncertainties refer to static α and θ. These errors are propagated from the error bars of the power spectrum.

The power law ranges from scales of 0.3 kpc (beam FWHM) to about 8 kpc. The fall-off in power at small scales is caused by instrumental beam effects. However, the larger scales are not greatly affected by the beam. To measure the spectral index we ignore the first spectral bin, which is possibly affected by map boundary limits as well as the last two beam affected bins. The dotted line in Fig. 8 represents a Kolmogorov slope of ϵ(k) ∝ k^-1.67 compared with the plot using θ = 45°. We achieve a perfect fit for the central part of the retrieved power spectrum, which is also inside the larger statistical variance of the parts at smaller k. There is a noticeable dependence of the spectral index on α with minimal and maximal slopes of 1.5 and 2.1, which should be regarded as improbable worst cases. θ has a big impact on the resulting magnetic field strengths. Assuming a larger θ implies a shorter interaction path with the magnetic field and therefore greater necessary magnetic field strengths. The derived B₀ in the cluster centre assumes that the magnetic profile given by Eq. (17) can be extrapolated with a constant α into the cluster centre. However, the sensitivity starts at regions from about 30 kpc from the cluster centre, see Fig. 4, and a flattening of the magnetic profile at smaller radii is well possible. B₅₀ gives the directly probed magnetic field strengths at 50 kpc distance.

The impact of the instrumental beam on the calculated magnetic field strength is low. A continuation of the power law to the small scales would give an additional 0.2 μG at most.

5. Comparison with other methods

The magnetic power spectrum of Hydra A was previously measured by Vogt & Enßlin (2005) and by Laing et al. (2008). Both find significantly lower central magnetic field strengths than we do. This is because of a much flatter magnetic field profile inferred characterised by a smaller index α. Moreover Vogt & Enßlin (2005) found a bump, a cut-off at spatial scales of about 2 kpc and more, but the spectral index for small scales is still similar to our current results. However, Laing finds a much flatter spectral index of about 0.8. Table 3 compares the results for a projection angle of θ = 45°.

To partly reconstruct the results of Vogt & Enßlin (2005), we applied the two simplifications they have used in the calculation of the covariance matrix in Eq. (18),

• 1.

  extending the second lower integration bound to –∞;

• 2.

  merging the two window functions f assuming the magnetic autocorrelation changes on much smaller scales than the window function.

Then we compared the results with the Fourier method of Vogt, but adding noise effects, the cavity structure, and a more advanced triangular parameterisation of the power spectrum, see Fig. 9. We reconstructed the spectral bump found by Vogt & Enßlin (2005) and found a perfect match between both methods. Figure 10 shows the difference to our method without these simplifications.

Consistent with Vogt & Enßlin (2005), the most likely α is shifted to much lower values if one adopts these simplifications, as can be seen in Fig. 11. However, instead of a falling plateau of the posterior between α = 0.1 and α = 0.8 and a cut-off for the unphysical range α < 0.1 we again find a Gaussian distribution, but with a maximum at 0.1. This very low α value indicates that the bump found by Vogt & Enßlin (2005) may be explained by a power leak on large scales produced by an implicit increase of the large-scale contribution of the window function.

One of the major simplifications used by Laing et al. (2008) is a linear and position-independent relation between the RM power spectrum and the magnetic autocorrelation⁴, which is only valid for a deep and statistically homogeneous Faraday screen. Both assumptions are not very well met for the Hydra A cluster. If we simulate a constant window function by setting α =  −1 (for which the effective window function, as given by Eq. (19), becomes flat) and use symmetric integration bounds in Eq. (18) we get a spectral index of about 1.1. For α = 0.25 we get a central magnetic field strength of 12.9 μG, which is much closer to their results.

The remaining difference to the spectral index of 0.8 and field strength of 19 μG may be partly caused by the beam correction applied by Laing et al. (2008), which lifts the high-frequency part of the spectrum, which is possibly contaminated with noise residuals from the RM-map-making algorithms. An implementation of a beam correction as a convolution in real space is not feasible in REALMAF, because it would be computationally too expensive and an implementation in Fourier space is in general only possible with simplifying assumptions. However, the exact frequency-dependent instrumental beam is unknown anyway, and we believe its effect can be neglected at scales larger than the beam size.

6. Conclusions

6.1. Summary

We developed a new method to retrieve magnetic power spectra from Faraday rotation maps implemented in the REALMAF (REAL MAgnetic Fields) code, where we modelled the magnetic autocorrelation function in real space.

Furthermore we introduced and verified a way to take the solenoidal character of the magnetic fields into account by using a simple full spherical model for the magnetic autocorrelation and multiplying with a proportionality factor. In contrast to the method by Vogt & Enßlin (2005) with a model in Fourier space, this permits us to alleviate some previously required simplifications in the processing. However, we are able to reconstruct their results by synthetically adding their simplifications to our method.

When we apply REALMAF to the north lobe of the Hydra A galaxy cluster, we find a power law power spectrum in the spatial range of 0.3 kpc to 8 kpc, with no visible large-scale turnover within this range, as reported by Vogt & Enßlin (2005). The magnetic density profile, which we modelled with an index α applied to the electron density profile, seems to follow the electron density profile with α = 1 ± 0.5. The somewhat large uncertainty on α also affects the spectral index of the power law, which is then in the range of 1.5 to 2.1. The exact projection angle of the used radio lobe Hydra A north is not well restricted, but it only seriously affects the absolute value of the magnetic energy density in the cluster and not the spectral slope. Depending on the angle we obtain field strengths ranging from 20 μG up to 80 μG in the cluster centre. However, these are extrapolations, because the measurement is only sensitive at a distance of about 50 kpc from cluster centre, where we find fields between 10 μG and 30 μG. The extrapolation of the model into the centre may be questionable.

Theories describing the generation of magnetic fields by turbulence like Subramanian et al. (2006), Enßlin & Vogt (2006) and Schekochihin & Cowley (2006) forecast a large-scale turnover, which might be visible inside our spectral range. Compared with Vogt & Enßlin (2005) we found no large-scale turnover within the spectral range up to 8 kpc. A turnover at larger scales is, however, not excluded by our analysis.

6.2. Outlook

The REALMAF code will permit us to analyse the Faraday rotation data of other galaxy clusters. This will hopefully allow us to study magnetic field spectra over a large range of spatial scales and in a variety of cluster environments. A usable map of Hydra A south would further increase the probed range at large scales and a map with a higher resolution would increase the range at small scales. We published the REALMAF code for general scientific usage under an open source license⁵.

Upcoming telescopes such as EVLA, LOFAR, LWA, ASKAP, and SKA will certainly provide the necessary high-fidelity Faraday rotation data of many more galaxy clusters for further studies with REALMAF. This may help to figure out the origin of cluster magnetic fields and the role they play in the metabolism of the intra-cluster medium.

---

Acknowledgments

Special thanks to Corina Vogt, who made her code used in Vogt & Enßlin (2005) available for insight, thanks to Daria Guidetti and Robert Laing for discussions and to an anonymous referee for constructive comments, and thanks to Henrik Junklewitz and Niels Oppermann for comments on the manuscript. This research was performed in the framework of the DFG Forschergruppe 1254 “Magnetisation of Interstellar and Intergalactic Media: The Prospects of Low-Frequency Radio Observations”.

References

Appendix A: Model functions

Here we describe the construction of our spectral model functions. We start with an ansatz for M_L(r) for the case $\mathrm{\nabla }\mathrm{·}\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{=}\mathrm{0}$ and use a Gaussian multiplied with a cosine. The corresponding M_N(r) is obtained by Eq. (23). Thus we adopt $\hspace{0.17em}\begin{array}{}\begin{array}{}\end{array}{\mathit{M}}_{\mathrm{L}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\frac{\mathrm{1}}{\sqrt{\mathrm{2}\mathit{\pi }}}{\mathit{b}}_{\mathit{i}}^{\mathrm{3}}\exp \mathrm{[}\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}\frac{{\mathit{r}}^{\mathrm{2}}}{{\mathit{a}}_{\mathit{i}}^{\mathrm{2}}}\mathrm{]}\cos \mathrm{\left[}{\mathit{b}}_{\mathit{i}}\mathit{r}\mathrm{\right]}\hspace{0.17em}\mathit{,}\hspace{0.17em}\mathrm{and}\begin{array}{}\end{array}\\ \begin{array}{}\end{array}{\mathit{M}}_{\mathrm{N}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\frac{{\mathit{b}}_{\mathit{i}}^{\mathrm{3}}}{\mathrm{2}{\mathit{a}}_{\mathit{i}}^{\mathrm{2}}\sqrt{\mathrm{2}\mathit{\pi }}}\exp \mathrm{\left[}\mathrm{-}\frac{{\mathit{r}}^{\mathrm{2}}}{\mathrm{2}{\mathit{a}}_{\mathit{i}}^{\mathrm{2}}}\mathrm{\right]}\mathrm{\left(}\mathrm{\right(}\mathrm{2}{\mathit{a}}_{\mathit{i}}^{\mathrm{2}}\mathrm{-}{\mathit{r}}^{\mathrm{2}}\mathrm{)}\cos \mathrm{\left[}{\mathit{b}}_{\mathit{i}}\mathit{r}\mathrm{\right]}\\ \begin{array}{}\end{array}\mathrm{+}{\mathit{a}}_{\mathit{i}}^{\mathrm{2}}{\mathit{b}}_{\mathit{i}}\mathit{r}\sin \mathrm{\left[}{\mathit{b}}_{\mathit{i}}\mathit{r}\mathrm{\right]}\mathrm{)}\mathit{.}\end{array}$(A.1)These functions are Fourier transformed analytically using

$\mathit{M}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\mathrm{4}\mathit{\pi }{\mathrm{\int }}_{\mathrm{0}}^{\mathrm{\infty }}{\mathit{r}}^{\mathrm{2}}\mathit{M}\mathrm{\left(}\mathit{r}\mathrm{\right)}\frac{\sin \mathrm{\left[}\mathit{kr}\mathrm{\right]}}{\mathit{kr}}\mathrm{d}\mathit{r,}$see Enßlin & Vogt (2003). We bind the free parameter a_i as

${\mathit{a}}_{\mathit{i}}\mathrm{=}\frac{\mathrm{3}}{\mathrm{4}}{\mathit{b}}_{\mathit{i}}\mathit{,}$because in this case ${\mathit{\omega }}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\mathrm{2}{\mathit{M}}_{\mathrm{N}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{+}{\mathit{M}}_{\mathrm{L}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{k}\mathrm{\right)}$ is fully positive and concentrated around a maximum at k = 3 b_i. This allows us to put the power spectrum together by relatively well separated and located frequency bands, as Fig. A.2 illustrates. Figure A.1 shows the model functions in position space. Especially M_N(r) has a strong negative part, because any closed field line has to come somewhere back through the plane perpendicular to it.

If we assume full isotropy, M_zz(r) becomes spherically symmetric. In this case a suitable model function can be obtained by Eq. (A.1): ${\mathit{M}}_{\mathit{zz}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\frac{\mathrm{1}}{\mathrm{3}}{\mathit{\omega }}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\frac{\mathrm{1}}{\mathrm{3}}\mathrm{\left(}{\mathit{M}}_{\mathrm{L}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\right(}\mathit{r}\mathrm{)}\mathrm{+}\mathrm{2}{\mathit{M}}_{\mathrm{N}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{(}\mathit{r}\mathrm{\left)}\mathrm{\right)}\mathit{.}$(A.2)This implicitly assumes M_N(r) = M_L(r). Thus we get

• a)

  using $\mathrm{\nabla }\mathrm{·}\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{=}\mathrm{0}$in position space: $\hspace{0.17em}\begin{array}{}\begin{array}{}\end{array}{\mathit{M}}_{\mathrm{L}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\frac{{\mathit{b}}_{\mathit{i}}^{\mathrm{3}}}{\sqrt{\mathrm{2}\mathit{\pi }}}\exp [\mathrm{-}\frac{\mathrm{8}}{\mathrm{9}}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}{\mathit{r}}^{\mathrm{2}}]\cos \mathrm{\left[}{\mathit{b}}_{\mathit{i}}\mathit{r}\mathrm{\right]}\begin{array}{}\end{array}\\ \begin{array}{}\end{array}{\mathit{M}}_{\mathrm{N}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\frac{{\mathit{b}}_{\mathit{i}}^{\mathrm{3}}}{\mathrm{18}\sqrt{\mathrm{2}\mathit{\pi }}}\exp [\mathrm{-}\frac{\mathrm{8}}{\mathrm{9}}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}{\mathit{r}}^{\mathrm{2}}]\mathrm{\left(}\mathrm{2}\mathrm{\right(}\mathrm{9}\mathrm{-}\mathrm{8}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}{\mathit{r}}^{\mathrm{2}}\mathrm{\left)}\cos \mathrm{\right[}{\mathit{b}}_{\mathit{i}}\mathit{r}\mathrm{]}\\ \begin{array}{}\end{array}\mathrm{-}\mathrm{9}{\mathit{b}}_{\mathit{i}}\mathit{r}\sin \mathrm{\left[}{\mathit{b}}_{\mathit{i}}\mathit{r}\mathrm{\right]}\mathrm{)}\mathit{,}\end{array}$(A.3)in Fourier space: $\hspace{0.17em}\begin{array}{}\begin{array}{}\end{array}{\mathit{M}}_{\mathrm{L}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\frac{\mathrm{27}}{\mathrm{64}}\mathit{\pi }(\frac{\mathit{k}\mathrm{-}{\mathit{b}}_{\mathit{i}}}{\mathit{k}}\exp [\mathrm{-}\frac{\mathrm{9}\mathrm{(}{\mathit{b}}_{\mathit{i}}\mathrm{-}\mathit{k}{\mathrm{)}}^{\mathrm{2}}}{\mathrm{32}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}}]\begin{array}{}\end{array}\\ \begin{array}{}\end{array}\mathrm{+}\frac{\mathit{k}\mathrm{+}{\mathit{b}}_{\mathit{i}}}{\mathit{k}}\exp [\mathrm{-}\frac{\mathrm{9}\mathrm{(}{\mathit{b}}_{\mathit{i}}\mathrm{+}\mathit{k}{\mathrm{)}}^{\mathrm{2}}}{\mathrm{32}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}}])\\ \begin{array}{}\end{array}{\mathit{M}}_{\mathrm{N}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\frac{\mathrm{27}}{\mathrm{128}}\mathit{\pi }((\frac{\mathrm{9}\mathrm{(}{\mathit{b}}_{\mathit{i}}\mathrm{-}\mathit{k}{\mathrm{)}}^{\mathrm{2}}}{\mathrm{16}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}}\mathrm{-}\mathrm{1})\exp [\mathrm{-}\frac{\mathrm{9}\mathrm{(}{\mathit{b}}_{\mathit{i}}\mathrm{-}\mathit{k}{\mathrm{)}}^{\mathrm{2}}}{\mathrm{32}{\mathit{b}}^{\mathrm{2}}}]\\ \begin{array}{}\end{array}\mathrm{+}(\frac{\mathrm{9}\mathrm{(}{\mathit{b}}_{\mathit{i}}\mathrm{+}\mathit{k}{\mathrm{)}}^{\mathrm{2}}}{\mathrm{16}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}}\mathrm{-}\mathrm{1})\exp [\mathrm{-}\frac{\mathrm{9}\mathrm{(}{\mathit{b}}_{\mathit{i}}\mathrm{+}\mathit{k}{\mathrm{)}}^{\mathrm{2}}}{\mathrm{32}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}}])\mathit{,}\end{array}$(A.4)and

• b)

  using full isotropy in position space: $\hspace{0.17em}\begin{array}{}\begin{array}{}\end{array}{\mathit{M}}_{\mathit{zz}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\frac{{\mathit{b}}_{\mathit{i}}^{\mathrm{3}}}{\mathrm{27}\sqrt{\mathrm{2}\mathit{\pi }}}\exp [\mathrm{-}\frac{\mathrm{8}}{\mathrm{9}}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}{\mathit{r}}^{\mathrm{2}}]\begin{array}{}\end{array}\\ \begin{array}{}\end{array}\mathrm{\times }\mathrm{(}\mathrm{(}\mathrm{27}\mathrm{-}\mathrm{16}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}{\mathit{r}}^{\mathrm{2}}\mathrm{)}\cos \mathrm{\left[}{\mathit{b}}_{\mathit{i}}\mathit{r}\mathrm{\right]}\mathrm{-}\mathrm{9}{\mathit{b}}_{\mathit{i}}\mathit{r}\sin \mathrm{[}{\mathit{b}}_{\mathit{i}}\mathit{r}{\mathrm{]}}^{\mathrm{)}}\mathit{,}\end{array}$(A.5)in Fourier space:

  Fig. A.2 Model functions ${\mathit{M}}_{\mathrm{N}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{k}\mathrm{\right)}$, ${\mathit{M}}_{\mathrm{L}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{k}\mathrm{\right)}$ and ω(i)(k) in Fourier space with ${\mathit{a}}_{\mathit{i}}\mathrm{=}\frac{\mathrm{3}}{\mathrm{4}}$ and bi = 1.

  $\hspace{0.17em}\begin{array}{}\begin{array}{}\end{array}{\mathit{M}}_{\mathit{zz}}^{\mathrm{\left(}\mathit{i}\mathrm{\right)}}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\frac{\mathrm{27}}{\mathrm{1024}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}\mathit{k}}\mathit{\pi }\exp [\mathrm{-}\frac{\mathrm{9}\mathrm{(}{\mathit{b}}_{\mathit{i}}\mathrm{+}\mathit{k}{\mathrm{)}}^{\mathrm{2}}}{\mathrm{32}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}}]\begin{array}{}\end{array}\\ \begin{array}{}\end{array}\mathrm{\times }\mathrm{(}\mathrm{16}{\mathit{b}}_{\mathit{i}}^{\mathrm{3}}\mathrm{+}\mathrm{9}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}\mathit{k}\mathrm{+}\mathrm{18}{\mathit{b}}_{\mathit{i}}{\mathit{k}}^{\mathrm{2}}\mathrm{+}\mathrm{9}{\mathit{k}}^{\mathrm{2}}\mathrm{+}\mathrm{9}{\mathit{k}}^{\mathrm{3}}\\ \begin{array}{}\end{array}\mathrm{+}\exp [\frac{\mathrm{9}\mathit{k}}{\mathrm{8}{\mathit{b}}_{\mathit{i}}}](\mathrm{-}\mathrm{16}{\mathit{b}}_{\mathit{i}}^{\mathrm{3}}\mathrm{+}\mathrm{9}{\mathit{b}}_{\mathit{i}}^{\mathrm{2}}\mathit{k}\mathrm{-}\mathrm{18}{\mathit{b}}_{\mathit{i}}{\mathit{k}}^{\mathrm{2}}\mathrm{+}\mathrm{9}{\mathit{k}}^{\mathrm{3}}))\mathit{.}\end{array}$(A.6)

Appendix B: Calculation of magnetic the field characteristics

Here we provide the used formulae for the magnetic field characteristics.

• a)

  Magnetic field strength The magnetic energy density can be calculated by integrating over the power spectrum or analogously taking the autocorrelation function in real space at the origin ω(0). The magnetic field strength in the centre of the cluster, where f = 1 and therefore $\begin{array}{}{}_{\mathrm{􏽥}}\\ \mathit{B}\end{array}\mathrm{=}\mathit{B}$, is then $\mathit{B}\mathrm{=}\sqrt{\mathit{\omega }\mathrm{\left(}\mathrm{0}\mathrm{\right)}}\mathit{.}$(B.1)

• b)

  Magnetic autocorrelation length ${\mathit{\lambda }}_{\mathit{B}}\mathrm{=}\frac{{\mathrm{\int }}_{\mathrm{-}\mathrm{\infty }}^{\mathrm{\infty }}\mathrm{d}\mathit{r\omega }\mathrm{\left(}\mathit{r}\mathrm{\right)}}{\mathit{\omega }\mathrm{\left(}\mathrm{0}\mathrm{\right)}}\mathrm{·}$(B.2)

• c)

  1-dimensional power spectrum Transformation from the three-dimensional power spectrum ω(k) to the one-dimensional ϵ(k) is done using (Enßlin & Vogt 2003) $\mathit{ϵ}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}\frac{{\mathit{k}}^{\mathrm{2}}\mathit{\omega }\mathrm{\left(}\mathit{k}\mathrm{\right)}}{\mathrm{2}\mathrm{(}\mathrm{2}\mathit{\pi }{\mathrm{)}}^{\mathrm{3}}}\mathrm{·}$(B.3)

• d)

  Spectral index using a Bayesian approach The spectral index can be calculated from the slope of the power spectrum in logarithmic scale. The fitting formula is: $\mathit{ϵ}\mathrm{\left(}\mathit{k}\mathrm{\right)}\mathrm{=}{\mathit{ϵ}}_{\mathrm{0}}(\frac{\mathit{k}}{{\mathit{k}}_{\mathrm{0}}}{)}^{\mathrm{-}\mathit{r}}\mathrm{=}\mathit{q}{\mathit{k}}^{\mathrm{-}\mathit{r}}\mathit{.}$The variance from the above power law is then defined as: $\mathit{\delta s}\mathrm{=}{\mathit{s}}_{\mathit{i}}\mathrm{-}\mathit{q}{\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}\mathit{.}$The allowed statistical deviations of s are imprinted in a covariance matrix D_ij =  ⟨ δs_i δs_j ⟩ , which is automatically obtained from the inverse of the Hessian matrix calculated with the maximum a posteriori approach as described in Sect. 3.3. The probabilities of the spectral amplitudes s of the power spectrum can be approximated as a Gaussian, where we can neglect the normalisation factor, because it is independent of q and r: $\mathit{P}\mathrm{\left(}\mathit{\delta s}\mathrm{\right|}\mathit{q,r}\mathrm{)}\mathrm{\propto }{\mathrm{e}}^{\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}}\mathit{\delta }{\mathit{s}}^{\mathit{T}}{\mathit{D}}^{-1}\mathit{\delta s}}\mathit{.}$Only a part of the power spectrum is used to fit the power law. D_ij is then the projection (cut) of the full inverse Hessian matrix $\begin{array}{}{\mathrm{􏽥}}_{}\\ {\mathit{D}}_{\mathit{ij}}\end{array}$ obtained from Eq. (16) ${\mathit{D}}_{\mathit{ij}}\mathrm{=}\mathrm{⟨}\mathit{\delta }{\mathit{s}}_{\mathit{i}}\mathit{\delta }{\mathit{s}}_{\mathit{j}}\mathrm{⟩}\mathrm{=}{\mathit{\psi }}^{\mathit{T}}\begin{array}{}\mathrm{􏽥}\\ {\mathit{D}}_{\mathit{ij}}\end{array}\mathit{\psi ,}$with the projection operator ψ =  { 0,...,1,...1,...,0 }  selecting the used spectral bins for the fit. Using the Bayesian theorem from Eq. (13): $\mathit{P}\mathrm{\left(}\mathit{q,r}\mathrm{\right|}\mathit{\delta s}\mathrm{)}\mathrm{\propto }\mathit{P}\mathrm{(}\mathit{\delta s}\mathrm{\left|}\mathit{q,r}\mathrm{\right)}\mathit{,}$where we already applied a uniform prior distribution P(r,q) = const. P(q,r | δs) has to be maximised to achieve the desired slope r within a maximum-a-posteriori approach. To numerically find the optimum, the gradients and the 2 × 2 Hessian matrix are necessary: $\begin{array}{ccc}\frac{\mathit{\partial }\ln \mathit{P}}{\mathit{\partial r}}& \mathrm{=}& \mathrm{-}\mathit{q}\ln \mathrm{\left[}{\mathit{k}}_{\mathit{i}}\mathrm{\right]}{\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}{\mathit{D}}_{\mathit{ij}}^{-1}\mathrm{(}{\mathit{s}}_{\mathit{j}}\mathrm{-}\mathit{q}{\mathit{k}}_{\mathit{j}}^{\mathrm{-}\mathit{r}}\mathrm{)}\\ \frac{\mathit{\partial }\ln \mathit{P}}{\mathit{\partial q}}& \mathrm{=}& {\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}{\mathit{H}}_{\mathit{ij}}^{-1}\mathrm{(}{\mathit{s}}_{\mathit{j}}\mathrm{-}\mathit{q}{\mathit{k}}_{\mathit{j}}^{\mathrm{-}\mathit{r}}\mathrm{)}\\ \frac{{\mathit{\partial }}^{\mathrm{2}}\ln \mathit{P}}{\mathit{\partial }{\mathit{q}}^{\mathrm{2}}}& \mathrm{=}& \mathrm{-}{\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}{\mathit{D}}_{\mathit{ij}}^{-1}{\mathit{k}}_{\mathit{j}}^{\mathrm{-}\mathit{r}}\\ \frac{{\mathit{\partial }}^{\mathrm{2}}\ln \mathit{P}}{\mathit{\partial }{\mathit{r}}^{\mathrm{2}}}& \mathrm{=}& \mathit{q}\mathrm{\left(}\mathit{q}\ln \mathrm{\right[}{\mathit{k}}_{\mathit{i}}\mathrm{\left]}{\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}{\mathit{D}}_{\mathit{ij}}^{-1}\ln \mathrm{\right[}{\mathit{k}}_{\mathit{j}}\mathrm{]}{\mathit{k}}_{\mathit{j}}^{\mathrm{-}\mathit{r}}\\ & & \mathrm{-}\mathrm{\left(}\ln \mathrm{\right[}{\mathit{k}}_{\mathit{i}}\mathrm{\left]}{\mathrm{)}}^{\mathrm{2}}{\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}{\mathit{D}}_{\mathit{ij}}^{-1}\mathrm{\right(}{\mathit{s}}_{\mathit{j}}\mathrm{-}\mathit{q}{\mathit{k}}_{\mathit{j}}^{\mathrm{-}\mathit{q}}\mathrm{\left)}\mathrm{\right)}\\ \frac{{\mathit{\partial }}^{\mathrm{2}}\ln \mathit{P}}{\mathit{\partial q\partial r}}& \mathrm{=}& \mathrm{2}\mathit{q}\ln \mathrm{\left[}{\mathit{k}}_{\mathit{i}}\mathrm{\right]}{\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}{\mathit{D}}_{\mathit{ij}}^{-1}{\mathit{k}}_{\mathit{j}}^{\mathrm{-}\mathit{r}}\mathrm{-}\ln \mathrm{\left[}{\mathit{k}}_{\mathit{i}}\mathrm{\right]}{\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}{\mathit{D}}_{\mathit{ij}}^{-1}{\mathit{s}}_{\mathit{j}}\mathit{.}\end{array}$An expression for the second unknown q can be found analytically. If the initial value of r is given, the corresponding initial value q can be calculated easily. This is necessary because the first guess must be close to the maximum to attain convergence. $\begin{array}{ccc}\mathit{q}\mathrm{=}\frac{{\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}{\mathit{D}}_{\mathit{ij}}^{-1}{\mathit{s}}_{\mathit{j}}}{{\mathit{k}}_{\mathit{i}}^{\mathrm{-}\mathit{r}}{\mathit{D}}_{\mathit{ij}}^{-1}{\mathit{k}}_{\mathit{j}}^{\mathrm{-}\mathit{r}}}\mathrm{·}& & \end{array}$

All Tables

All Figures