Issue 
A&A
Volume 563, March 2014



Article Number  A56  
Number of page(s)  15  
Section  Numerical methods and codes  
DOI  https://doi.org/10.1051/00046361/201322825  
Published online  06 March 2014 
Detecting and quantifying stellar magnetic fields
Sparse Stokes profile approximation using orthogonal matching pursuit
LeibnizInstitut für Astrophysik (AIP),
An der Sternwarte 16,
14482
Potsdam,
Germany
email: tcarroll@aip.de; kstrassmeier@aip.de
Received: 10 October 2013
Accepted: 15 November 2013
Context. In recent years, we have seen a rapidly growing number of stellar magnetic field detections for various types of stars. Many of these magnetic fields are estimated from spectropolarimetric observations (Stokes V) by using the socalled centerofgravity (COG) method. Unfortunately, the accuracy of this method rapidly deteriorates with increasing noise and thus calls for a more robust procedure that combines signal detection and field estimation.
Aims. We introduce an estimation method that provides not only the effective or mean longitudinal magnetic field from an observed Stokes V profile but also uses the net absolute polarization of the profile to obtain an estimate of the apparent (i.e., velocity resolved) absolute longitudinal magnetic field.
Methods. By combining the COG method with an orthogonalmatchingpursuit (OMP) approach, we were able to decompose observed Stokes profiles with an overcomplete dictionary of waveletbasis functions to reliably reconstruct the observed Stokes profiles in the presence of noise. The elementary wave functions of the sparse reconstruction process were utilized to estimate the effective longitudinal magnetic field and the apparent absolute longitudinal magnetic field. A multiresolution analysis complements the OMP algorithm to provide a robust detection and estimation method.
Results. An extensive MonteCarlo simulation confirms the reliability and accuracy of the magnetic OMP approach where a mean error of under 2% is found. Its full potential is obtained for heavily noisecorrupted Stokes profiles with signaltonoise variance ratios down to unity. In this case a conventional COG method yields a mean error for the effective longitudinal magnetic field of up to 50%, whereas the OMP method gives a maximum error of 18%. It is, moreover, shown that even in the case of very small residual noise on a level between 10^{3} and 10^{5}, a regime reached by current multiline reconstruction techniques, the conventional COG method incorrectly interprets a large portion of the residual noise as a magnetic field, with values of up to 100 G. The magnetic OMP method, on the other hand, remains largely unaffected by the noise, regardless of the noise level the maximum error is no greater than 0.7 G.
Key words: stars: magnetic field / Sun: magnetic fields / radiative transfer / line: formation / line: profiles / polarization
© ESO, 2014
1. Introduction
There are various ways of obtaining information about the presence of magnetic fields in stellar atmospheres like Ca II H & K or Xray measurements. But the most direct way is certainly provided via the Zeeman effect, i.e., the magnetically induced line splitting and the associated line polarization. This has opened the way for essentially two different approaches to measuring the magnetic fields of cool stars. Unpolarized Zeeman broadening measurements from Stokes I profiles has led to a large amount of information about the magnetic fields of cool stars over the past decades (e.g., Robinson 1980; Saar 1988; Valenti et al. 1995; JohnsKrull & Valenti 1996; JohnsKrull 2007). Because the Zeeman splitting (and broadening) is directly proportional to the magnetic field strength, Zeeman broadening measurements provide good estimates of the underlying average surface field. In the stellar astrophysical literature, the term magnetic flux for the product of surface filling factor and magnetic field strength (f·B) is often used, but we prefer to denote it more accurately as the mean or average magnetic field strength, since f·B is neither the (relative) magnetic flux through the stellar surface nor a measure of the magnetic flux density. Polarization or spectropolarimetric measurements, on the other hand, provide additional information about the orientation of the field vector and therefore allow, in principle, measurement of real flux densities.
There are at least two obstacles that may prevent the full applicability of spectropolarimetric measurements to cool stars. First, the magnetic fields are relatively small in size and magnitude, which leads to very weak circular polarization (Stokes V) and even lower linear polarization (Stokes Q and U) signals. Second, the contributions from different polarities tend to mutually cancel out the circular polarization signal. A circumstance that mitigates the situation to some extent is rapid rotation, which partly lifts the mutual cancellation of the Stokes V spectra. This has led to the socalled ZeemanDoppler imaging (ZDI) or MagneticDoppler imaging (MDI) methods (Semel 1989; Donati et al. 1997; Piskunov & Kochukhov 2002), which allows reconstructing, in an inversion approach, the entire surface distribution of the magnetic field vector from phaseresolved spectropolarimetric observations. This approach has contributed significantly to our understanding of cool stars magnetic fields and stellar activity in general in recent years (e.g., Donati et al. 2003, 2006, 2008; Kochukhov et al. 2004; Petit et al. 2008; Hussain et al. 2009; Morin et al. 2010; Carroll et al. 2012). However, what is good for the ZDI approach, namely rapid rotation, has an adverse effect on the Zeeman broadening approach and vice versa, such that both approaches are complementary to some degree; see also Reiners (2012) for a comprehensive overview of measuring stellar magnetic fields on cool stars.
The present paper focuses on spectropolarimetric observations that do not deliver a series of phaseresolved spectra but rather single snapshotlike observations from either rapidly or slowly rotating active stars. This is the typical situation for stars with long rotation periods and/or in situations where the amount of telescope time is limited. Such observations are responsible for a number of exciting magnetic field detection in recent times (e.g., Aurière et al. 2007, 2009; Lignières et al. 2009; Aurière et al. 2010; Grunhut et al. 2010; KonstantinovaAntova et al. 2010; Sennhauser & Berdyugina 2011; Petit et al. 2011; Fossati et al. 2013). Because many of the observed polarimetric line profiles are buried in noise, the investigations rely on a prior line profile reconstruction technique, such as the leastsquaredeconvolution (LSD; Donati et al. 1997), principalcomponentanalysis (PCA; Carroll & Kopf 2007), or singularvaluedecomposition (SVD; Carroll et al. 2012). However, even after these preprocessing steps, the extracted line profiles exhibit a fraction of noise that makes their interpretation not always straightforward in terms of a reliable detection and magnetic field estimation.
Inspired by the work of Asensio Ramos & Lopez Ariste (2010) who propose a compressive sensing framework for spectropolarimetry, we utilize the compressibility of polarimetric signals to approximate the observed Stokes V profiles by a sparse decomposition with a wavelet frame. This sparse linear representation facilitates a simple and noiserobust magnetic detection and estimation method^{1}. For rapidly rotating stars, the sparse approximation of the observed signal with basis functions allows not only the effective mean longitudinal magnetic field to be measured but also the absolute value of the resolved (i.e. apparent) longitudinal magnetic field.
The paper is organized as follows. In Sect. 2 we describe the weakfield approximation and derive the centerofgravity method for diskintegrated observations to define the effective longitudinal magnetic field. Furthermore, we also define the apparent longitudinal magnetic field that describes the diskintegrated absolute value of the longitudinal magnetic field component as measured from the net absolute circular polarization of the Stokes V profile. In Sect. 3 we introduce the concept of a sparse Stokes profile approximation and magnetic field detection. We describe in some detail the orthogonal matching pursuit (OMP) algorithm and the signal dictionary, which provides the building blocks for the Stokes profile approximation. At the end of this section, we combine the approximation method with a detection algorithm to obtain reliable estimates for the magnetic quantities. Section 4 gives an illustration of the magnetic OMP method and highlights the benefit of the complementary definition of the effective and apparent longitudinal magnetic field. In Sect. 5 we investigate the sparsity of Stokes profiles and demonstrate that with the described dictionary a sparse approximation of the Stokes profile is possible for a wide parameter range. An extensive numerical assessment and evaluation of the accuracy of the magnetic OMP method are also performed in Sect. 5, along with a comparison with the conventional centerofgravity method. A summary is presented in Sect. 6.
2. Magnetic field estimation
Our starting point will be the socalled centerofgravity (COG) method (Semel 1967; Rees & Semel 1979). Before we describe the COG method, we briefly introduce the concept of the weakfield approximation from which we later derive the effective and apparent longitudinal magnetic field.
2.1. DiskIntegrated weakfield approximation
The local weakfieldapproximation (WFA; e.g., Landi degl’Innocenti 1992; Stenflo 1994) makes the assumption that the atmosphere of a resolved region is permeated by a homogeneous, weak, and depthindependent magnetic field. Then, a proportionality exists between the observed local Stokes V profile and the spectral derivative of the observed local Stokes I profile, which can be derived from the polarized transfer equation using a firstorder Taylor expansion of the line profile function (e.g., Jefferies et al. 1989). This proportionality or linear relation can be written as (1)where γ is the angle between magnetic field vector and the lineofsight (LOS), g the effective Landé factor, and λ_{B} the Zeeman splitting expressed in wavelength units is given by (2)For the sake of brevity, we denote by λ the relative wavelength shift from the line center λ_{0} hereafter and assume that the Stokes profiles are normalized by the local continuum intensity. We moreover introduce the variable μ = cosγ and to write Eq. (1) as (3)where B_{l} denotes the longitudinal component of the magnetic field vector. The proportionality between the observed local Stokes V and the observed local Stokes I profile depends on the wavelength. In the local case, i.e. resolved (solar) observations, this dependency is usually neglected since the same scaling factor, i.e. αBμ, applies for all wavelengths. This condition is not generally true for unresolved stellar observations.
Given that the assumptions for the weakfield approximation are applicable, it is readily seen that for the local and the resolved case, one can obtain an estimate for the longitudinal component B_{l} of the magnetic field through (4)For each spectral line, a good estimate up to which magnetic field strength the WFA remains a valid approximation, can be deduced from the socalled Zeeman saturation (Stenflo 1994). This is the regime from which on the Stokes V amplitudes no longer grow linearly but instead the individual sigma components begin to shift apart according to Eq. (2). The Zeeman saturation regime is reached for the majority of Zeeman sensitive spectral lines in the optical at around one kiloGauss.
In stellar spectropolarimetric observations, the only observables are the diskintegrated Stokes profiles. The obvious requirement for applying the weakfield approximation to stellar spectra is the existence of a wavelength independent linear scaling between the diskintegrated Stokes V profile, denoted hereafter as , and the spectral derivative of the diskintegrated Stokes . For the following disk integration, we choose a Cartesian coordinate system where the zaxis points along the LOS, the yaxis is in the plane defined by the LOS and the rotation axis, and the xaxis is perpendicular to that plane. Using Δ λ_{Rot} = (λ_{0} v sin i)/c, the diskintegrated weakfield approximation for a rotating star can then be written in its general form as (5)where we introduced the spatialdependent longitudinal magnetic field strength B_{l}(x,y) = B(x,y)μ(x,y).
The observed Stokes vector depends on the magnetic field, whereas the local intensity has no explicit dependency on the magnetic field under the weakfield assumption. The intensity profile originates over the entire observable disk. This disconnection, in general, prevents applying the WFA to unresolved stellar observation because the diskintegrated profile shape of the spectral derivative of Stokes does not coincide with the observed Stokes profile. However, the assumption that the thermodynamic parameters in the atmosphere and across the observable disk do not change between magnetic and fieldfree regions, as well as the assumption that the longitudinal magnetic field on the stellar surface has no explicit spatial dependence, makes the WFA applicable also to diskintegrated spectra. In this case one can obtain a similar relation to the one in Eq. (4) by defining a mean longitudinal field strength ⟨B_{l}⟩ , which allows one to separate the magnetic field from the integral in Eq. (5). The diskintegrated WFA can then compactly written as (6)This linear scaling relation must hold for all wavelengths λ. If this is not the case, owing to a spatial dependence of the underlying magnetic field and rapid rotation, the WFA can no longer be applied. The applicability of the diskintegrated WFA is also limited to fields below the Zeeman saturation regime. Because the disk integration is a linear operation, the Zeeman saturation regime, in general, is determined not by the mean longitudinal field but by the local field strengths on the surface of the star.
2.2. The effective longitudinal magnetic field
The centerofgravity method does in fact bypass the limitation that the observable Stokes profile must be proportional to the spectral derivative of the Stokes profile by using the firstorder moment of the Stokes profile (Mathys 1989). This firstorder moment of the Stokes V signal can be obtained by using Eq. (5), (7)where Λ denotes the integration limits that cover the entire line profile. Instead of deriving the COG method from the weakline limit we use the WFA as a starting point (Semel 1967). Under the WFA, the profile of a spectral line does not alter its shape, and we may introduce a linedependent intensity distribution function η, which accounts for geometrically induced radiative transfer effects (e.g., limb darkening and other possible temperature effects). We then write, for the Stokes I profile at position x,y on the stellar disk, (8)Substituting this into Eq. (7) gives (9)By taking advantage of the Cartesian coordinate system, where values along the xaxis experience different degrees of Doppler shifts Δλ, we can perform the following change of variables, x = Δλ/Δλ_{Rot}, to write (10)We define the total fluxweighted longitudinal magnetic field distribution along a small a strip dx or d(Δλ) of equal velocity as (11)This allows us to write Eq. (10) as a convolution integral (12)or in a more compact way (13)The diskintegrated Stokes profile is also given by the convolution (14)where η(Δλ) is defined by (15)With this definition we may separate a mean longitudinal magnetic field from Eq. (11) that we call the effective longitudinal magnetic field B_{eff}. Using Eqs. (14) and (15), we can rearrange Eq. (12) to obtain the following expression for the effective longitudinal magnetic field (16)Performing the integration by parts in the denominator and keeping in mind that we implicitly deal with normalized profiles, we obtain (17)where W_{I} is the equivalent width of the diskintegrated Stokes I profile. This is the COG method for diskintegrated observations introduced for solar observations by Semel (1967); Rees & Semel (1979).
Making a transformation of the relative wavelength according to λ → vλ_{0}/c, we obtain the following equation for the velocity domain (18)where the central wavelength λ_{0} is given in Å. This form is extensively used in spectropolarimetric observations where the spectral line profiles are often preprocessed by transforming them into the velocity domain during a reconstruction process (see e.g., Donati et al. 1997; Wade et al. 2000).
2.3. The apparent longitudinal magnetic field
If we now consider the case where the projected rotational velocity of the star is sufficiently large and the effective longitudinal magnetic field is zero owing to a balance of positive and negative polarities, the diskintegrated Stokes V signal of the star does not necessarily cancel out and vanishes. In fact, the net absolute circular polarization; i.e., the integral of the absolute value of the observed Stokes V signal over wavelength is in general not zero. Even though we measure no effective field, the polarized profile is clearly detectable and apparently a magnetic field must be present. Of course, this is one of the effects that is utilized by ZDI to resolve surface magnetic fields on rapidly rotating stars. Having such a nonvanishing net absolute circular polarization means that the magnetic field exhibits a flux imbalance for each or some of the resolved surface areas on the stellar disk, even though the overall diskaveraged magnetic field is flux balanced.
To characterize this imbalance, we define the apparent longitudinal magnetic field. We start by defining the integrated Stokes V^{∗} along a small isoradial velocity strip d(Δλ) as (19)The diskintegrated firstorder moment of the Stokes V profile can then be written as (20)The total amount of the effective longitudinal magnetic field can be expressed by (21)This is, of course, not an observable, but assuming that we have a number of n discrete and spectroscopically resolvable radialvelocity strips across the surface of the star that produce separated Stokes V profiles, we may write (22)where is the Stokes V profile originating in the ith resolved radialvelocity strip. In this case the observed Stokes V profile is the composition of n individual contributions coming from different velocityresolved regions. Using Eqs. (22) and (20), we can finally define the apparent longitudinal magnetic field as (23)or in velocity coordinates (24)The apparent longitudinal magnetic field measures the maximum resolved field that can be obtained from the observation. In general, we have , and the amount of the B_{app} depends on the intrinsic line width of the local Stokes I profile. However, it can give an estimate for the distribution of the absolute effective field over the stellar disk. Although it is not clear in the first place how to obtain the apparent longitudinal field and the flux density distribution from the diskintegrated Stokes V profile, we present an easy and fast way for estimating this quantity in the next section.
3. Sparse approximation of Stokes profiles and magnetic field detection
The basic idea of our proposed approach is to find an accurate approximation of the observed Stokes V profile by decomposing the original signal with a set of well chosen elementary waveforms. The stagewise approximation by these elementary basis functions allows us to apply the equations for the effective and apparent longitudinal magnetic field to each of these elementary functions and eventually obtain the desired magnetic quantities. The decomposition we describe in the following is based on a prescribed overcomplete set of elementary functions, called signal or profile atoms. Overcomplete here means that the number of elements provided by the dictionary is redundant, i.e. larger than the actual dimension of the signal space. In contrast to many nonredundant (e.g., orthogonal) transformations, a linear expansion of a given signal profile with an overcomplete set of signal atoms often facilitates a more efficient and sparse approximation (Tropp 2004; Donoho et al. 2006). By using a suitable and redundant set of profile atoms, we show in the following that a sparse approximation of observed Stokes line profiles, allows a resolved analysis of the effective and apparent longitudinal magnetic flux density.
3.1. Orthogonal matching pursuit
The matching pursuit (MP) algorithm is an iterative algorithm for adaptive signal reconstruction and approximation. A given signal or line profile is decomposed into a linear expansion of dictionary elements. The actual realization of the dictionary is discussed in Sect. 3.2. For the moment it suffices to consider the dictionary as a collection of possibly lineardependent elementary waveforms. We begin with a description of the original MP algorithm introduced by Mallat & Zhang (1993) and then of its improvement called OMP (Pati et al. 1993), which we use for our analysis.
Following Mallat & Zhang (1993) let us define a dictionary D(x) as a collection of a large number of elementary wave functions [x_{1},x_{2},...,x_{n}] formally defined in Hilbert space ℋ where each vector is of unit norm ∥ x_{i} ∥ = 1. A given signal f ∈ ℋ is approximated by MP in a first step by projecting f onto a vector x_{l} ∈ D such that (25)where Rf is the current residual of the approximation. Since the residual Rf of the current estimate is orthogonal to x_{l}, we have the following energy conservation (26)To minimize ∥ Rf ∥ ^{2} the vector x_{l} ∈ D is chosen such that the maximum projection (i.e., correlation) is found between the residual vector and all dictionary atoms, (27)The algorithm now iteratively decomposes the residual Rf by repeated projections onto the dictionary atoms. The recursive algorithm can be written in a compact way if we set the initial values for the residual R^{0}f = f and the initial approximation to f_{0} = 0, (28)An increasingly closer approximation to the signal is obtained by repeating steps (I) to (III). A proof of convergence for the MP algorithm can be found in Mallat & Zhang (1993).
The MP algorithm can be improved by realizing that the steps taken by the MP are not optimal in the sense that a newly chosen dictionary vector is not orthogonal to the previously selected vectors (Pati et al. 1993). This can be resolved by orthogonalizing the newly selected vector relative to the n−1 previously selected vectors. The orthogonalization, which is essentially a GramSchmidt procedure (Bronstein et al. 1997), uses an auxiliary vector z that expresses the yet unexplained component of the newly selected vector x_{l}. For the nth iteration we may write (29)where the initial condition is set to z_{0} = x_{l0}. The new residual R^{n + 1}f can then be expressed by the projection of the current residual R^{n}f on z_{n}, (30)using the relation in Eq. (29). The latter can also be written as (31)Compared to the ordinary MP algorithm a component is now subtracted in a direction that is orthogonal to all previously selected vectors. Beginning with the initial conditions R^{0}f = f and f_{0} = 0 as well as z_{0} = x_{l0} the recursive algorithm for the OMP can then be expressed as (32)After n iterations we obtain the following sparse approximation of our signal, (33)A convergence analysis and more results for the OMP algorithm can be found in Davis et al. (1997). The set of selected signal atoms { x_{l0},x_{l1},...,x_{ln} }, its corresponding projection weights, and auxiliary vectors give us the opportunity to analyze any given signal in terms of a small number of dictionary basis functions.
3.2. The wavelet dictionary for Stokes profile approximation
To facilitate a sparse and compact representation of a given Stokes profile, we need a set of basis functions that closely resemble the individual building blocks of a diskintegrated circular polarized profile. These building blocks are the local Stokes V profiles originating in a small resolved surface area of the star. From the definition of a general spectral line profile (Gray 2005), we know that a Gaussian function can provide a reasonably good approximation of a spectral absorption profile in cases where damping is not too strong. Under the regime of the WFA, where the Stokes V profile is proportional to the derivative of the intensity profile, it seems obvious that an appropriate and easy way to approximate a composite (i.e., diskintegrated) Stokes V profile can be achieved by using a dictionary of the derivatives of Gaussian functions. In fact, we use a wavelet frame constructed from scaled and translated versions of first derivatives of a Gaussian (Torrence & Compo 1998). We define the dictionary D as a set of elementary functions or atoms Ψ_{j,k} of unit norm as (34)The atoms are constructed by a set of scaled and translated versions of a mother wavelet Ψ. For the mother wavelet we choose, as mentioned above, the first derivative of a Gaussian function, (35)The wavelet function obtained by a scaled and translated version of the mother wavelet can be expressed in the wavelength domain as (36)where s_{j} and λ_{k} denote the scaling and translation parameter. To obtain a discrete set of wavelet functions, we use the following discretization of the scale parameter (37)where s_{0} defines the smallest resolvable scale, which we set to 2δλ, twice the wavelength step of the observations or synthetic calculations. The parameter Δ_{r} is a variable value that provides the resolution of the transform and which is set throughout this paper to Δ_{r} = 0.125. This provides a reasonable tradeoff between resolution and computational cost (e.g., Torrence & Compo 1998). The largest scale L is determined according to (38)where N is the number of wavelength points in the line profile. With the wavelet function Eq. (36), we can express the wavelet coefficients of an observed spectrum V(λ) as the dot product between the Stokes V profile and the wavelet function, (39)where the integration limits R are large enough to cover the entire profile. The projections (i.e., dot products) of Eq. (33) used in the OMP algorithm can be computed easily by calculating the convolution between the current residual and the wavelet function for all scales s_{j}.
The decomposition by a restricted set of wavelets with different scales provides a multiscale representation of our original signal profile. We use this multiscale decomposition to apply an approach called multiresolution support (Starck & Murtagh 2006) to complement the OMP algorithm with an efficient detection procedure.
3.3. Magnetic field detection
We first start with the assumption that our observation vector (i.e., spectral line profile) V contains the true signal vector S and additive noise N such that we may write, for the observed signal, (40)Separating the signal from noise can be cast in the framework of nonparametric signal estimation. One very successful way of reconstructing an unknown signal from noisy observation is thresholding introduced and theoretically investigated, e.g., by Donoho (1995). For our detection analysis we take a similar thresholding approach that uses the multiscale decomposition of the signal to test the significance of every wavelet coefficient individually. The basic strategy of shrinkage or thresholding methods is to consider the observed signal in a transformed (e.g., wavelet) domain rather than in the original data domain. Similar to Fourier filtering, waveletbased detections methods rest on the idea that the noise contribution exhibits a different statistical behavior in the transformed domain.
Fig. 1
Orthographic maps of the random magnetic distribution of the test star at four different rotational phases φ. 

Open with DEXTER 
The OMP algorithm approximates the observed line profile in an iterative process by projecting the dictionary atoms onto the current residual of the observation while picking the most correlated signal atom in each iterative cycle. Using a wavelet dictionary, the projections can be efficiently calculated by a convolution, i.e., wavelet transform of the observed line profile. These projections between an individual wavelet function of scale j at position k and the current residual of the observation V in the presence of noise can be written thanks to the linearity of the wavelet transform as (41)which we may write according to Eq. (37) as the wavelet coefficient w_{j,k}(42)When including the position index k into a vector notation, this can be written in a compact way by using an operator or matrix notation (43)where W_{V},W_{S}, and W_{N} represent the convolution over the wavelength index k on a specific scale j. The expectation value for the wavelet transform of the noise contribution is E { W_{N} } = 0, and the covariance is given by , where T denotes the transpose. Uncorrelated white noise reduces the covariance matrix Σ to a diagonal matrix. Unless the matrix is orthogonal (i.e., by using an orthogonal set of basis functions), we obtain for each scale j a different noise level σ_{j}. Before we address the problem of estimating the noise level on each resolution scale j, we take a closer look at the detection problem. The wavelet transform yields a level or resolutiondependent representation of the observed spectral line profile in terms of the wavelet coefficients w_{j,k}. A decision about whether this wavelet coefficient is significant, i.e., whether it includes signal information or not, can be cast in a hypothesistesting framework. For this, we can state the null hypothesis H_{0} such that w_{j,k} contains only noise and no signal information. Whether the null hypothesis will be rejected or not depends on probability P_{n}, (44)where τ is a detection threshold. The null hypothesis will be rejected if P_{n} is smaller than a given significance level ϵ; i.e. P_{n}(τ) < ϵ . For a Gaussian noise distribution with zero mean and standard deviation σ, we may write the probability density for w_{j,k} as (45)The probability P_{n} of rejecting H_{0} can then be calculated by integrating over all weights w_{j,k} that are greater than a threshold τ, (46)Since the wavelet transform of the Stokes profile can have both positive and negative values, we also need to integrate over all negative weights w_{j,k} up to the threshold value −τ(47)When assuming stationary noise and choosing a specific significance level ϵ, this reduces to the following decision or threshold detection rule (e.g., Starck & Murtagh 2006) (48)where k depends on the value of ϵ. Choosing ϵ = 0.0027 results in k = 3, the well known 3σ detection threshold. This thresholding detection scheme can be introduced readily into the OMP algorithm of Sect. 3.1 by implementing the threshold detection rule of Eq. (48) as an additional step (Ib) after calculating the maximum projection (i.e. inner product) in step (I). If the currently selected wavelet coefficient is not significant, it will be deleted from the dictionary for this iteration cycle, and the algorithm restarts the current iteration cycle with step (I). The noise level σ_{j} for each scale j can be obtained by a wavelet transform of simulated noise. This process yields leveldependent thresholds and has the advantage that different types of noise (e.g., Gaussian or Poisson), as well as correlated noise, can be modeled.
4. Application to synthetic observations
We now use the OMP algorithm to estimate the effective and apparent longitudinal magnetic field. Before we begin with a statistical analysis of the accuracy of our proposed method, we give an illustrative example to highlight the functionality of the algorithm. The decomposition of an observed Stokes V profile into the building blocks of the wavelet dictionary provides the opportunity for each selected signal or profile atom to be individually analyzed in terms of the effective and apparent longitudinal magnetic field. Thanks to the linearity of the algorithm, we can evaluate Eqs. (17) and (23) for each signal atom and iteratively approximate the desired magnetic quantities. Each contributing signal atom found by the OMP algorithm identifies a resolved structure in the wavelength (or velocity) domain. The overall number of iterations then also determines the summation index n for the apparent longitudinal field in Eq. (23). The stopping criterion is given by the detection threshold as soon as the approximation enters the noise and no more significant signal atoms are detected.
In the following example we have built a stellar surface model with our iMap code (Carroll et al. 2012). We created a random distribution of a radial magnetic field; i.e., each surface segment has a Gaussian random distribution of its radial magnetic field strength with a mean value centered at zero Gauss and a standard deviation of 500 Gauss. The effective temperature of the test star is 5250 K, and it has solar abundance and a gravity of log (g) = 4.0. The projected rotational velocity vsini is 35 km s^{1}, and micro and macroturbulence were set to 2 km s^{1}. The inclination of the rotational axis is 90°. Because the original surfacegrid resolution with a 5° by 5° segmentation creates a strong cancellation of the Stokes V signal we decided, for this example case, to additionally smear out the distribution over the surface with a Gaussian filter that has an angular spread of 15°.
Fig. 2
Stokes V of the Zeeman sensitive iron line FeI 6173 for the random field distribution of the test star at phase 0.0. 

Open with DEXTER 
Fig. 3
Magnetic OMP approximation at different iteration cycles. The true synthetic profile (solid line) and the approximation of the OMP method (dashed line). 

Open with DEXTER 
This smoothing leads to large random clusters of magnetic fields over the surface. The resulting magnetic field surface distribution is shown in Fig. 1. For the synthetic calculation, we used the magnetically sensitive iron line at 6173 Å. Line parameters were taken from the VALD line database (Piskunov et al. 1995; Kupka et al. 1999). The resulting synthetic Stokes V profile for phase 0.0 is shown in Fig. 2. To highlight the performance of the method, we used noisefree Stokes profiles in this first test. The stopping criterion for the noisefree case was substituted by a convergence criterion where the iteration was stopped as soon as there was no significant improvement in the rootmeansquared (RMS) error of the approximation. The approximations are illustrated in Fig. 3 where four snapshots at different stages in the approximation process are shown. The magnetic OMP algorithm gradually identifies the most coherent signal atoms on different scales and positions before it finally approximates the entire observed Stokes V profile with a linear combination of all selected dictionary atoms. As can be seen in Fig. 3, the algorithm finds the one bestmatching signal atom in the first iteration (upper left). In the approximations 5 (upper right), 10 (lower left), and 20 (lower right) one can follow the rapid improvement after adding more and more signal atoms to the residuals. Already at iteration cycle 10 (i.e., 10 dictionary atoms), the differences between the original synthetic observation and the approximation is hardly visible. In Fig. 4 the RMS error of the approximation is plotted over the iteration cycle, which also demonstrates the rapid convergence of the magnetic OMP algorithm. However, the approximation of the pure Stokes V profile is not the main task of the magnetic OMP method here, because it is supposed to give accurate estimates of the effective and apparent magnetic field as well.
Fig. 4
The root mean square approximation error over the iteration cycle. 

Open with DEXTER 
The effective and apparent longitudinal magnetic field is estimated during the approximation process from each contributing signal atom. To compare the estimated values with the true values, we need to extract the true effective longitudinal magnetic field from the model star. This is straightforward and only requires projecting the magnetic field vector of each surface segment onto the lineofsight and eventually to sum the areaweighted contribution up from all visible surface segments. The extraction of the true apparent longitudinal magnetic field value from the model deserves some more explanation. In principle we could use the same process as for the effective longitudinal field and sum the absolute value up instead of the plain values from the stellar model, but this would not take the cancellation of the Stokes V signals coming from within isoradial velocity strips into account. We would therefore overestimate the apparent field; in fact, we would retrieve the total absolute effective longitudinal magnetic field. To obtain an adequate value for the apparent longitudinal magnetic field, we need to bring the surface distribution of the longitudinal magneticflux density into an appropriate isoradial velocity binning before we can add up their absolute contributions. The resolution of this isoradial velocity binning depends on the width of the used (local) spectral line; i.e., the broader the spectral line, the less surface flux can be deduced, and vice versa. Therefore we first translate the surface distribution of the longitudinal field of the model star onto a finebinned rotational velocity coordinate axis. The resolution of this onedimensional coordinate axis has the same resolution as our synthetic spectral line profile, i.e. 0.5 km s^{1}. The distribution of the longitudinal magnetic flux density along the rotation velocity is then convolved with an areanormalized local Stokes I profile to obtain the measurable longitudinal fluxdensity spectrogram over the rotational velocity which is depicted in Fig. 5.
We then sum over the absolute values of the longitudinal fluxdensity spectrogram to finally obtain the apparent longitudinal magnetic field that is compared to the estimates of the OMP method. In our test case (phase = 0.0), the true effective longitudinal magnetic field inferred from the model is 4.67 Gauss and the apparent longitudinal magnetic field is 16.43 Gauss. From the magnetic OMP algorithm, we obtain a value of 4.53 Gauss for the effective longitudinal magnetic field and 16.78 for the apparent longitudinal magnetic field. An exhaustive statistical evaluation is given in the next section, but one can already see that besides the good approximation and the rapid convergence, the method also yields accurate values for both magnetic quantities.
There is another interesting effect that highlights the complementary nature of the effective and apparent longitudinal magnetic field. At phase 0.25 of our test star, the polarities in the longitudinal magnetic field are almost perfectly balanced over the visible hemisphere. The resulting profile (again synthesized for the iron line FeI 6173) is shown in Fig. 6. The almost perfect balancing of polarities is not obvious from the line profile itself. However, if we estimate the longitudinal magnetic field by the COG method’s Eq. (16) and the magnetic OMP algorithm, we obtain 0.005 G or 0.007 G, respectively, for the effective longitudinal magnetic field.
Fig. 5
Distribution of the longitudinal magnetic flux density across the rotational velocity for phase 0.0 of the test star. 

Open with DEXTER 
Fig. 6
Stokes V profile of the Zeeman sensitive iron line FeI 6173 for the random distribution at phase 0.25 of the test star. 

Open with DEXTER 
We thus have an extremely small longitudinal magnetic field or none at all. This is not what one would expect from the mere visual inspection of the Stokes V profile in Fig. 6, which apparently shows a clear net absolute polarization. In fact, the amplitude of the Stokes V signal is even greater than the one of phase = 0.0, shown in Fig. 2. This demonstrates the strength of the definition of the apparent longitudinal magnetic field; the magnetic OMP algorithm detects a clear apparent longitudinal field of 18.63 G (true value 17.98 G), which is slightly stronger even than in phase 0.0. The quantity of the apparent magnetic longitudinal field is therefore of particular interest in cases where balanced field distributions cause the firstorder moment of the Stokes V profile to vanish.
Parameter ranges used for the model atmosphere and synthetic line calculations.
5. Numerical experiments and statistical evaluation
In this section we assess the accuracy of the OMP algorithm under a broad range of model conditions and for different spectral lines. For that reason we synthesized a large number of Stokes I and Stokes V profiles for different magnetic surface distributions and atmospheric conditions. The model parameters are the effective temperature, metallicity, surface gravity, projected rotational velocity, and the surface magnetic field and its distribution. For each set of the atmospheric parameters, we created a random surface distribution in the same way as for the previous example in Fig. 1. The atmospheric parameters were randomly chosen from within an interval given in Table 1. For the model atmospheres, we chose Kurucz/Atlas9 models (Castelli & Kurucz 2004). Because the model atmospheres are provided on a fixed grid of step sizes of 250 K for the effective temperature, 0.5 for the logarithmic abundance and gravity, the atmospheres are interpolated for each randomly chosen set of atmospheric parameters. As Zeemansensitive spectral lines we chose three magnetically sensitive lines; the iron lines FeI 5497 (g_{eff} = 2.22), FeI 6173 (g_{eff} = 2.50), and FeI 8468 (g_{eff} = 2.50). All line parameters were again taken from the VALD line database (Piskunov et al. 1995; Kupka et al. 1999). For each spectral line, we created a sample of 5000 randomly chosen atmospheric parameters, rotational velocities, and magnetic surface distributions to eventually synthesize a set of corresponding Stokes I and Stokes V profiles with the forward module of our iMap code.
In total, we calculated a set of 15 000 Stokes I and Stokes V profile, which were then analyzed by our magnetic OMP algorithm. The true effective longitudinal magnetic field is directly extracted from the magnetic surface distribution of the synthetic model star. The true apparent longitudinal magnetic field is obtained from the convolved magnetic spectrograms as described above. Because the individual trials, sets of atmospheric parameters, and surface magnetic field values have different scales, we chose a relative performance metric to describe the accuracy.
The overall performance is judged by the mean absolute percentage error, while for illustrating the error distribution we use the relative percentage error. The relative percentage error is defined as (49)where is the true magnetic field, and the calculated magnetic field of the ith sample. Finally, the mean absolute percentage error (MAPE) is defined by (50)
5.1. Sparsity of Stokes profiles
Before we evaluate the performance of the OMP algorithm, we look more closely at the general sparsity of Stokes profiles. A sparse representation is often achieved by transforming a vector from the original data domain into a domain where the transformed coefficients (i.e., expansion coefficients) allow a more compact representation of the information (e.g., Fourier or wavelet transform). A vector is called sparse if there is a representation where most of the vector entries are zero or close to zero. A sparse approximation is performed by restricting the sparse transformation to the largest expansion coefficients such that there is no great loss of signal information.
Recently, Asensio Ramos & Lopez Ariste (2010) have shown that Stokes profiles are sparse and compressible under various transformations (e.g., wavelets and empirical basis functions). Compressible here means that the values of the sorted expansion coefficients exhibit an exponential decay that in turn results in a small approximation error (Candes & Wakin 2008). Although we saw from the numerical simulation of the last section that a sparse representation of a synthetic Stokes profile is possible with our dictionary elements, we also want to test empirically if this also holds for the entire set of our 15 000 Stokes profiles. To quantify the approximation error, we use the relative error, where is the approximation and S^{∗} the original Stokes vector.
Fig. 7
The mean relative approximation error for the OMP dictionary and the PCA expansion. The relative error of the approximation is plotted over the number of signal atom (eigenprofiles for the PCA). The error declines steeply for both expansion methods; however, the OMP expansion shows a better performance, i.e., sparsity behavior compared to the PCA. 

Open with DEXTER 
Fig. 8
Error distribution of the longitudinal magnetic field calculated from 15 000 synthetic Stokes profiles by the magnetic OMP method (left) and the centerofgravity (COG) method (right). Both methods show very good accuracy with a mean absolute percentage error of 1.87% (magnetic OMP) and 1.79% (COG). 

Open with DEXTER 
To compare the performance of the OMP expansion relative to a known sparse decomposition, we also calculated a principal component analysis (PCA) of the entire synthetic data set. A PCA decomposes a set of observed Stokes profiles into an empirical set of orthogonal eigenprofiles (e.g., Carroll et al. 2008). As has been shown by Asensio Ramos & Lopez Ariste (2010) and Asensio Ramos & Allende Prieto (2010), a PCA expansion can result in a compact and sparse representation of Stokes profiles. Figure 7 shows how the dictionary used with the OMP algorithm performs against a PCA decomposition of our synthetic data set. Both curves (solid OMP, dashed PCA) show a rapid decline of the mean approximation error with increasing numbers of signal atoms or eigenvectors, respectively. The performance of our dictionary of Gaussian derivatives is better than that of the empirical basis functions (i.e., eigenprofiles) obtained by the PCA. The reason for that is twofold: First, the signal atoms of our dictionary are selected a priori to match the building blocks (i.e., elementary waveforms) of our problem and second the redundancy of our dictionary facilitates a more efficient expansion then that of the orthogonal eigenvectors of the PCA.
In contrast to the empirical basis formed by the PCA, the OMP dictionary is not required to be a set of orthogonal basis functions, and this generally provides better flexibility and adaptivity to approximate a given signal (Candes & Wakin 2008).
To reduce the relative approximation error with the OMP algorithm below 5%, it takes on average only nine signal atoms. Using 22 signal atoms the relative error is on average smaller than 1%, which demonstrates that our dictionary allows a sparse representation of Stokes profiles for the parameter regime given in Table 1.
5.2. The noisefree case
To obtain an idea about the principal accuracy of the magnetic OMP algorithm and its performance relative to the conventional COG method, we began with a noisefree test case where the entire test sample of 15 000 Stokes I, and Stokes V profiles were used without any noise contribution.
In this simulation the overall error value for the effective longitudinal magnetic field obtained from the magnetic OMP algorithm yielded a MAPE of 1.87%. For comparison the effective longitudinal magnetic field calculated by the COG method shows a similar MAPE of 1.79%. The same error for the apparent longitudinal magnetic field is 4.67%. The error for the three individual spectral lines show no apparent deviation from the overall error.
For the subset of synthetic Stokes V profiles of the FeI 5497, we obtain a MAPE for the effective (apparent) longitudinal magnetic field of 1.80% (4.57%), for the FeI 6173 lines a MAPE of 1.93% (4.69%), and for the FeI 8468 lines a MAPE of 1.88% (4.72%). The error distribution for the relative percentage error of the effective longitudinal field calculated by the OMP and the COG method is shown in Fig. 8. Both the OMP method and the conventional COG method provide equally good results. The error distribution of the apparent longitudinal magnetic field is shown in Fig. 9. Given that the apparent field is even harder to determine due to possible cancellation effects of the Stokes V signals, the error is still remarkably low, and the OMP method is able to recover the apparent field from the Stokes V profiles with good accuracy.
5.3. The noise case
Fig. 9
Error distribution of the apparent longitudinal magnetic field calculated by the magnetic OMP method. The mean absolute percentage error from the entire set of 15 000 synthetic test calculations is 4.67%. 

Open with DEXTER 
Fig. 10
Mean absolute percentage error (MAPE) of the effective longitudinal magnetic field calculated over the relative noise level for the COG method (open circles) and the magnetic OMP method (solid circles). 

Open with DEXTER 
Fig. 11
Mean absolute percentage error (MAPE) of the apparent longitudinal magnetic field calculated versus the relative noise level. 

Open with DEXTER 
In the following we test the accuracy of the magnetic OMP method and the COG method for the more relevant case when the Stokes V profiles are contaminated with noise. For that reason we have added an increasing amount of white noise to the 15,000 Stokes V profiles. Since the Stokes V profiles generated from the random surface distribution are of different magnitudes, we define a relative noise level η that is given by the relative variance of the Stokes V profile and that of the noise according to (51)where σ(S) and σ(N) are the standard deviations of the Stokes V profile and the noise contribution, respectively. Figure 10 shows how the relative error increases with higher noise levels. It is interesting how rapidly the accuracy of the COG method deteriorates (upper curve) compared to the relative robustness of the magnetic OMP method (lower curve). Already for a noise level of η_{rel} = 0.3, the accuracy of the COG method is significantly affected, whereas the OMP method still provides good accuracy with errors of less than 10%. Even for a noise level of unity, the OMP method has an error of just 18% where the COG method with an error of almost 50% can no longer be considered a meaningful tool for quantifying the magnetic field. For the apparent longitudinal magnetic field, the OMP method shows a very robust performance as well. The error distribution over the relative noise level is illustrated in Fig. 11. The good performance and robustness of the magnetic OMP method compared to the COG method can be understood by considering the noisefree Stokes V profile in Fig. 11. Because the OMP algorithm gradually approximates the noisy profile with a linear combination of dictionary profile atoms, it calculates the effective and apparent longitudinal magnetic field for each of the selected profile atoms separately. The detection mechanisms of Sect. 3.3 prevent the algorithm from approximating insignificant profile features and thus largely avoids the quantification of the contributing noise. Figure 12 shows one of the sample profiles generated for a random surface distribution of the magnetic field. The noiseless profile is approximated well by the OMP method, and the estimated effective and apparent longitudinal magnetic field is approximated with –41.04 G and 59.88 G, respectively, which is very close to the true values of –40.11 G and 61.65 G. The noiseless profile is approximated by the maximum number of 40 iterations, which means that 40 dictionary profile atoms are used to approximate the example profile. However, already the first three profile atoms are able to approximate the profile to such a degree that 95% of the signal variance are described. In Fig. 13 the same profile can be seen with varying degrees of noise contributions overplotted again by the approximation from the OMP algorithm. Despite the increasing relative noise level, the OMP algorithm finds a good approximation to the original Stokes V profile. The iteration for a noise level of η_{rel} = 0.25 stops already at seven iterations (i.e., 7 profile atoms). For η_{rel} = 0.5, 0.75, and 1.0, the iteration is stopped at cycles 6, 4, and 3, respectively. At a relative noise level of unity only three signal atoms are used to approximate the Stokes V profile. However, as one can see from the noiseless case, these three signal atoms already capture the essential profile shape and amplitude of the original Stokes V signal. Since the magnetic OMP algorithm interprets and analyzes each signal atom separately, the good approximation directly translates to a good estimates of the underlying effective and apparent longitudinal magnetic field. In this particular case, the estimations for the effective (apparent) longitudinal field from the magnetic OMP algorithm is –42.54 G (64.27 G) for a noise level of 0.25, –44.26 G (65.71 G) for a noise level 0.5, –36.04 G (68.85 G) for a noise level 0.75, and –34.77 G (72.54 G) for noise level of 1.0.
Fig. 12
A noiseless (noise level 0.0) example profile (solid line), and its OMP approximation (dashed line). 

Open with DEXTER 
Fig. 13
Magnetic OMP approximations (dashed lines) to a noisy Stokes V profile (solid lines) for different relative noise levels. Despite the increasing relative noise, the OMP approximation identifies the main characteristic and shape of the underlying noisefree profile. 

Open with DEXTER 
5.4. Absolute limits of the COG and OMP method
Until now we have analyzed the performance of the magnetic OMP and COG method by using relative noise levels. To shed more light on the absolute performance limits of the two methods, we take a closer look on the absolute values of the noise. Multiline reconstruction techniques like LSD (Donati et al. 1997) or SVD (Carroll et al. 2012) allow us to drastically lower the noise levels of spectropolarimetric observations. The thus increased signaltonoise (S/N) levels allows detecting faint signal profiles that were otherwise deeply buried in the noise.
A natural question that arises in the context of magnetic field estimation is how much of the residual noise is falsely interpreted as magnetic field and up to which limits we can reliably estimate a true underlying magnetic field. To address this problem, we could simply perform a firstorder perturbation to the COG method and apply standard error propagation to estimate the influence of the noise. However, for weak magnetic fields, the noise contribution relative to the true signal amplitude is in general not small anymore. We therefore simulate the process with a large and statistically significant number of synthetic profiles to determine the impact of typical absolute noise levels on the magnetic field estimation.
We assume additive Gaussian noise with zero mean. One might expect that this would lead to a vanishing value in the field estimation by the COG method. This is, however, not true for a finite spectral resolution. The random fluctuations and the unequal weighting in the velocity (or wavelength) domain can lead to relatively large spurious contributions to the firstorder moment, hence to the field estimation as shown in the following simulation. Note that any form of correlated noise would even exacerbate the problem since it would introduce a systematic bias in the evaluation of the firstorder moment.
With our assumption of white noise, a Stokes profile can be written as the sum of the true noiseless Stokes profile and a noise vector N, (52)The effective or mean longitudinal field as measured by the COG method can then be written as (53)We assume that the noise level is small enough that the estimation of the Stokes I equivalent width remains largely unaffected.
The firstorder moment of a noisy Stokes profile, as well as the magnetic field estimation, can then be split into a signal () and a noise () contributing part, (54)Written in this form we may ask to which noise level the following relation is valid, , or at which point does a noiseinduced magnetic field interfere with the magnetic field estimation from the true signal profile.
For the simulations we defined the absolute noise level η_{abs} as the standard noise deviation σ_{N} relative to the continuum of the normalized intensity profile, i.e., η_{abs} = σ_{N}, which is also commonly used to describe the S/N level (1/σ_{N}) in stellar polarimetry (Donati et al. 1997). To mimic a multiline, reconstructed line profile (e.g., LSD profile), we created a fictitious iron line with a mean Landé factor of 1.2 at a rest wavelength of 5000 Å. The oscillator strength and van der Waals damping were adjusted to give an equivalent width for the Stokes I profile of 100 mÅ. The Stokes I and Stokes V noise profiles were calculated with a resolution of λ/Δλ = 100 000, which is common to current highresolution observations. The spectral range covers a velocity of plus/minus 50 km s^{1} around the line center, and integration limits were set accordingly.
Fig. 14
The effective magnetic field response over the absolute noise level (logarithmic scale) for the COG method (open circles) and the OMP method (solid circles). The response curve shows how much of the noise is erroneously attributed to a magnetic field. On the left side, the response is shown for a noise range between 10^{5} and 10^{3}. The noise response of the COG method rapidly increases with the noise level from 2 to 100 G, whereas the OMP method remains largely unaffected by the noise with a flat response not larger than 0.7 G. On the right side, a closeup of the magnetic response curve for the region between 10^{5} to 10^{4} is shown. Even in such a low noise regime, the COG method falsely interprets a significant portion of the noise as magnetic field. The OMP method stays essentially flat at a level close to zero over the entire region. The gray shaded area marks the region of incorrect response for the COG method, i.e. , and defines the lower limit for magnetic field detections 

Open with DEXTER 
We chose noise levels from 10^{3} down to 10^{5} that correspond to a polarimetric sensitivity of 0.001%, a level that is reached by current reconstruction methods (Donati & Landstreet 2009). We selected 20 noise levels between 10^{5} and 10^{3}, equally spaced on a logarithmic scale. For each absolute noise level, we synthesized 10 000 noisy line profiles such that we have a total set of 200 000 simulated profiles. Owing to the linear superposition of the signal and noise contribution in Eq. (54), we do not depend on an underlying magnetic field and only need to simulate, besides the Stokes I profile, the noise part of the profile. To quantify the magnetic response to the noise, we calculated the effective longitudinal magnetic field for each synthetic noise profile with the COG and the OMP methods. To obtain a statistical meaningful average response for each of 20 noise levels, we computed the mean of the absolute effective field from each of the corresponding 10 000 field values.
In Fig. 14, we plotted the magnetic response of the COG and of the OMP method over the various noise levels. From the response curve, we can immediately identify how much of the noise is falsely interpreted as a magnetic field. As the noise level increases, the COG method is increasingly affected by the noise, and more and more of the content of each noise profile is interpreted as an effective longitudinal magnetic field. For the COG method, a noise level of 10^{3} results, on average, in an effective longitudinal magnetic field of 101 G. On the right side of Fig. 14, we expand the region between 10^{5} to 10^{4}. A few times 10^{5} is the region where many of the recently detected weak magnetic fields are reported (e.g., Aurière et al. 2009, 2010; Lignières et al. 2009; Grunhut et al. 2010; KonstantinovaAntova et al. 2010; Petit et al. 2011; Fossati et al. 2013). Despite the low noise level, we see that the COG method interprets a nonnegligible amount of the noise as a magnetic field (between 2 and 10 G).
The shaded areas in Fig. 14 mark the region of magnetic field values that can no longer be properly interpreted by the COG method, i.e., the region where . In that sense the response curves set the lower limit or a noise threshold for the estimation of weak magnetic fields. However, even when the true field strength is above these threshold values, the relative contribution from the noise can severely compromise the field estimation of the COG method as was shown in Sect. 5.3. For the magnetic OMP method, on the other hand, we see from Fig. 14 that the response is largely flat.
The magnetic field estimation from the OMP method is almost unaffected by the noise and shows only a small increase with the noise level. The incorrectly attributed field values are no higher than 0.7 G. This has its origin in the multiresolution thresholding scheme implemented in the OMP algorithm (see Sect. 3.3).
Only those signal or noise features that reach a certain significance are interpreted and quantified. This again results in the good performance and robustness against the contributing noise. A dark gray area that shows the incorrect response region of the OMP method is not visible in Fig. 14 due to the low noise response of the OMP method. The specific values for the noise response of the COG method also depend of course on the line parameters, such as wavelength, equivalent width, and Landé factor. However, the deviations from our results are expected to be small because the spread in the derived average parameters for multiline reconstruction techniques that use line lists of many hundreds or thousands of lines are relatively small. The spectral resolution is also a contributing factor to the error of the COG method since the statistical fluctuations that enter into the evaluation of the firstorder moment directly depend on the number of available wavelength points, such that lower spectral resolutions on average lead to larger errors (i.e., to a stronger noise response). However, a simulation and critical assessment of all these contributing factors in the COG estimation, are not subjects of the current paper.
6. Summary
This work presents a novel technique for detecting and quantifying stellar magnetic fields. By employing a sparse profile approximation in terms of an orthogonal matching pursuit algorithm, we were able to provide a robust method for detecting and estimating the effective longitudinal magnetic field.
By introducing the concept of the apparent longitudinal magnetic field, which describes the maximum of the resolvable absolute longitudinal surface field over the stellar disk, we provided a complementary measure to the effective or mean longitudinal magnetic field. For rapidly rotating stars, the definition of the apparent longitudinal magnetic field allows a quantification of the field even in situations where a smallscale and balanced surface field would otherwise result in a vanishing mean longitudinal magnetic field.
It was shown by an extensive numerical simulation with our iMap code (Carroll et al. 2012) that the accuracy of the OMP method is remarkably good with a mean absolute percentage error of only 1.87% for the effective longitudinal magnetic field and 4.67% for the apparent longitudinal magnetic field. However, the real strength of the new approach is its robustness against noise. Even for relative noise levels of unity, i.e., when the variance of the noise has the same magnitude as the actual Stokes V signal, the OMP method has a mean error of only 18%, whereas the conventional estimation by the COG method yields a mean error of almost 50%.
In an effort to understand the limitations of the conventional COG method for estimating weak magnetic fields from faint and noisy Stokes V signals, we simulated noise profiles in a lownoise regime from an absolute noise level of 10^{3} down to 10^{5}. This is a regime that is currently reached by multiline reconstruction techniques such as LSD (Donati et al. 1997) or SVD (Carroll et al. 2012) and will be reached by the next generation of spectropolarimeters like PEPSI at the 2 × 8.4 m Large Binocular Telescope (LBT) (Strassmeier et al. 2008). In these lownoise regimes, however, a nonnegligible fraction of the noise was incorrectly interpreted as a magnetic field by the COG method. Depending on the noise level, these falsely identified field values are between 2 and 100 G, making any attempt to estimate weak fields of a similar strength extremely difficult. The magnetic OMP method, on the other hand, again shows a remarkable robustness over the entire range of noise levels, where only 0.1 to 0.7 G are incorrectly attributed to a magnetic field. This makes the magnetic OMP method the method of choice for estimating weak magnetic fields from noisy Stokes profiles.
From a numerical perspective, the algorithm can be easily implemented and the iterative process is fast to evaluate^{2}. The magnetic OMP algorithm provides a viable tool for detecting and quantifying magnetic fields from spectropolarimetric observations. Although the focus in this work has been placed on diskintegrated Stokes profiles of stellar magnetic fields, this method is also directly applicable to resolved solar magnetic field observations.
An IDLCode of the magnetic OMP method is available under http://www.aip.de/People/tcarroll and http://www.aip.de/Members/tcarroll
We provide an IDLCode of the magnetic OMP method under http://www.aip.de/People/tcarroll/ and http://www.aip.de/Members/tcarroll/
Acknowledgments
The authors would like to thank the anonymous referee for the constructive and helpful comments, which helped to improve this manuscript.
References
 Asensio Ramos, A., & López Ariste, A. 2010, A&A, 509, A49 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Asensio Ramos, A., & Allen de Prieto, C. 2010, ApJ, 719, 1759 [NASA ADS] [CrossRef] [Google Scholar]
 Auer, L. H., &Heasley, J. N. 1978, A&A, 64, 67 [NASA ADS] [Google Scholar]
 Aurière, M., Wade, G. A.,Silvester, J., et al. 2007, A&A, 475, 1053 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Aurière, M., Wade, G. A.,KonstantinovaAntova, R., et al. 2009, A&A, 504, 231 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Aurière, M., Donati, J.F., KonstantinovaAntova, R., et al. 2010, A&A, 516, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Bronstein, I. N., Semendjajew, K. A., Musiol, G., & Mühlig, H. 1997, Taschenbuch der Mathematik (Frankfurt am Main: Verlag Harri Deutsch) [Google Scholar]
 Candes, E. J., &Wakin, M. B. 2008, IEEE Signal Process. Mag., 25, 21 [NASA ADS] [CrossRef] [Google Scholar]
 Carroll, T. A., &Kopf, M. 2007, A&A, 468, 323 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Carroll, T. A., Kopf, M., Ilyin, I., &Strassmeier, K. G. 2007, AN, 328, 1043 [NASA ADS] [CrossRef] [Google Scholar]
 Carroll, T. A., Kopf, M., &Strassmeier, K. G. 2008, A&A, 488, 781 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Carroll, T. A.,Strassmeier, K. G., Rice, J. B., & Künstler, A. 2012, A&A, 548, A95 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Castelli, F., & Kurucz, R. L. 2004, IAU Symp. 210, poster A20 [Google Scholar]
 Davis, G., Mallat, S., &Avellaneda, M. 1997, J. Constr. Approx., 13, 57 [CrossRef] [Google Scholar]
 Donati, J.F., &Landstreet, J. D. 2009, Annu. Rev. Astron. Astrophys., 47, 333 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Donati, J.F., Semel, M., &Praderie, F. 1989, A&A, 225, 467 [NASA ADS] [Google Scholar]
 Donati, J.F., Semel, M., &Rees, D. E. 1992, A&A, 265, 669 [NASA ADS] [Google Scholar]
 Donati, J.F., Semel, M., Carter, B. D., Rees, D. E., & CollierCameron, A. 1997, MNRAS, 291, 658 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Donati, J.F., Collier Cameron, A.,Hussain, G. A. J., & Semel, M. 1999, MNRAS, 302, 437 [NASA ADS] [CrossRef] [Google Scholar]
 Donati, J.F., Cameron, A., & Collier, Semel, M. et al. 2003, MNRAS, 345, 1145 [NASA ADS] [CrossRef] [Google Scholar]
 Donati, J.F., Howarth, I. D.,Jardine, M. M., et al. 2006, MNRAS, 370, 629 [NASA ADS] [CrossRef] [Google Scholar]
 Donati, J.F., Moutou, C., Fares, R., et al. 2008, MNRAS, 385, 1179 [NASA ADS] [CrossRef] [Google Scholar]
 Donoho, D. L. 1995, IEEE Trans. Inform. Theory, 41, 613 [CrossRef] [MathSciNet] [Google Scholar]
 Donoho, D. L., Tsaig, J., Drori, I., &Starck, J.L. 2006, IEEE Trans. Inform. Theory, 58, 1094 [CrossRef] [Google Scholar]
 Fossati, L., Kochukhov, O., Jenkins, J. S., et al. 2013, A&A, 551, A85 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gray, D. F. 2005, The Observation and Analysis of Stellar Photospheres, 3rd edn. (Cambridge University Press) [Google Scholar]
 Grunhut, J. H., Wade, G. A., Hanes, D. A., &Alecian, E. 2010, MNRAS, 408, 2290 [NASA ADS] [CrossRef] [Google Scholar]
 Hubrig, S., Schöller, M., Kholtygin, A. F., et al. 2012, A&A, 546, L6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hussain, G. A. J., Collier Cameron, A.,Jardine, M. M., et al. 2009, MNRAS, 398, 189 [NASA ADS] [CrossRef] [Google Scholar]
 Jefferies, J., Lites, B. W., &Skumanich, A. 1989, ApJ, 343, 920 [NASA ADS] [CrossRef] [Google Scholar]
 JohnsKrull, C. M. 2007, ApJ, 664, 975 [NASA ADS] [CrossRef] [Google Scholar]
 JohnsKrull, C. M., &Valenti, J. A. 1996, ApJ, 459, L95 [NASA ADS] [CrossRef] [Google Scholar]
 Kochukhov, O., Bagnulo, S., Wade, G. A., et al. 2004, A&A, 414, 613 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 KonstantinovaAntova, R., Aurière, M., Charbonnel, C., et al. 2010, A&A, 524, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kupka, F., Piskunov, N., Ryabchikova, T. A.,Stempels, H. C., & Weiss, W. W. 1999, A&AS, 138, 119 [NASA ADS] [CrossRef] [EDP Sciences] [MathSciNet] [PubMed] [Google Scholar]
 Landi degl’Innocenti, E. 1992, Solar Observations: Techniques and Interpretation (Cambridge University Press), 71 [Google Scholar]
 Landstreet, J. D. 1982, ApJ, 258, 639 [NASA ADS] [CrossRef] [Google Scholar]
 Landstreet, J. D., Silaj, J., Andretta, V., et al. 2008, A&A, 481, 465 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lignières, F., Petit, P., Böhm, T., &Aurière, M. 2009, A&A, 500, L41 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Martínez González, M. J., Asensio Ramos, A.,Carroll, T. A., et al. 2008, A&A, 486, 637 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mallat, S., &Zhang, S. 1993, IEEE Trans. Signal Process., 41, 3397 [NASA ADS] [CrossRef] [Google Scholar]
 Mathys, G. 1989, Fund. Cosm. Phys., 13, 143 [NASA ADS] [Google Scholar]
 Mathys, G. 2002, in Astrophysical Spectropolarimetry (Cambridge University Press), 101 [Google Scholar]
 Morin, J., Donati, J.F., Petit, P., et al. 2010, MNRAS, 407, 2269 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Pati, Y. C., Rezaiifar, R., & Krishnaprasa, P. S. 1993, Orthogonal Matching Pursuit: recursive function approximation with application to wavelet decomposition, Proc. of the 27th Annual Asilomar Conf. on Signals, Systems, and Computers, IEEE Comput. Soc. Press, 40 [Google Scholar]
 Petit, P., Dintrans, B., Solanki, S. K., et al. 2008, MNRAS, 388, 80 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Petit, P., Lignières, F., Aurière, M., et al. 2011, A&A, 532, L13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Piskunov, N., &Kochukhov, O. 2002, A&A, 381, 736 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Piskunov, N. E., Kupka, F., Ryabchikova, T. A., Weiss, W. W., &Jeffery, C. S. 1995, A&AS, 112, 525 [NASA ADS] [Google Scholar]
 Rees, D. E., &Semel, M. D. 1979, A&A, 74, 1 [NASA ADS] [Google Scholar]
 Reiners, A. 2012, Liv. Rev. Sol. Phys., 9, 1 [Google Scholar]
 Robinson, R. D., Jr. 1980, ApJ, 239, 961 [NASA ADS] [CrossRef] [Google Scholar]
 Saar, S. H. 1988, ApJ, 324, 441 [NASA ADS] [CrossRef] [Google Scholar]
 Semel, M. 1967, Ann. Astrophys., 30, 513 [NASA ADS] [Google Scholar]
 Semel, M. 1989, A&A, 225, 456 [NASA ADS] [Google Scholar]
 Sennhauser, C., &Berdyugina, S. V. 2011, A&A, 529, A100 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Skumanich, A., & López Ariste, A. 2002, ApJ, 570, 379 [NASA ADS] [CrossRef] [Google Scholar]
 Starck, J.L., & Murtagh, F. 2006, Astronomical Image and Data Analysis, Astron. Astrophys. Lib., 2nd edn. (SpringerVerlag) [Google Scholar]
 Starck, J.L., Murtagh, F., & Bijaoui, A. 1998, Image Processing and Data Analysis: The Multiscale Approach (Cambridge: University Press) [Google Scholar]
 Stenflo, J. O. 1994, Solar Magnetic Fields: Polarized Radiation Diagnostics, Astrophys. Space Sci. Lib. (Dordrecht; Boston: Kluwer Academic Publishers), 189 [Google Scholar]
 Strassmeier, K. G., Woche, M., Ilyin, I., et al. 2008, Proc. SPIE, 7014 [Google Scholar]
 Torrence, C., &Compo, G. P. 1998, Bull. Am. Meteor. Soc., 79, 61 [NASA ADS] [CrossRef] [Google Scholar]
 Tropp, J. A. 2004, IEEE Trans. Inform. Theory, 50, 2231 [CrossRef] [MathSciNet] [Google Scholar]
 Tropp, J. A. 2007, IEEE Trans. Inform. Theory, 53, 4655 [CrossRef] [Google Scholar]
 Valenti, J. A., Marcy, G. W., &Basri, G. 1995, ApJ, 439, 939 [NASA ADS] [CrossRef] [Google Scholar]
 Wade, G. A., Donati, J.F., Landstreet, J. D., &Shorlin, S. L. S. 2000, MNRAS, 313, 851 [NASA ADS] [CrossRef] [Google Scholar]
 Wade, G. A.,Auriére, M.,Bagnulo, S., et al. 2006, A&A, 451, 293 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
All Tables
All Figures
Fig. 1
Orthographic maps of the random magnetic distribution of the test star at four different rotational phases φ. 

Open with DEXTER  
In the text 
Fig. 2
Stokes V of the Zeeman sensitive iron line FeI 6173 for the random field distribution of the test star at phase 0.0. 

Open with DEXTER  
In the text 
Fig. 3
Magnetic OMP approximation at different iteration cycles. The true synthetic profile (solid line) and the approximation of the OMP method (dashed line). 

Open with DEXTER  
In the text 
Fig. 4
The root mean square approximation error over the iteration cycle. 

Open with DEXTER  
In the text 
Fig. 5
Distribution of the longitudinal magnetic flux density across the rotational velocity for phase 0.0 of the test star. 

Open with DEXTER  
In the text 
Fig. 6
Stokes V profile of the Zeeman sensitive iron line FeI 6173 for the random distribution at phase 0.25 of the test star. 

Open with DEXTER  
In the text 
Fig. 7
The mean relative approximation error for the OMP dictionary and the PCA expansion. The relative error of the approximation is plotted over the number of signal atom (eigenprofiles for the PCA). The error declines steeply for both expansion methods; however, the OMP expansion shows a better performance, i.e., sparsity behavior compared to the PCA. 

Open with DEXTER  
In the text 
Fig. 8
Error distribution of the longitudinal magnetic field calculated from 15 000 synthetic Stokes profiles by the magnetic OMP method (left) and the centerofgravity (COG) method (right). Both methods show very good accuracy with a mean absolute percentage error of 1.87% (magnetic OMP) and 1.79% (COG). 

Open with DEXTER  
In the text 
Fig. 9
Error distribution of the apparent longitudinal magnetic field calculated by the magnetic OMP method. The mean absolute percentage error from the entire set of 15 000 synthetic test calculations is 4.67%. 

Open with DEXTER  
In the text 
Fig. 10
Mean absolute percentage error (MAPE) of the effective longitudinal magnetic field calculated over the relative noise level for the COG method (open circles) and the magnetic OMP method (solid circles). 

Open with DEXTER  
In the text 
Fig. 11
Mean absolute percentage error (MAPE) of the apparent longitudinal magnetic field calculated versus the relative noise level. 

Open with DEXTER  
In the text 
Fig. 12
A noiseless (noise level 0.0) example profile (solid line), and its OMP approximation (dashed line). 

Open with DEXTER  
In the text 
Fig. 13
Magnetic OMP approximations (dashed lines) to a noisy Stokes V profile (solid lines) for different relative noise levels. Despite the increasing relative noise, the OMP approximation identifies the main characteristic and shape of the underlying noisefree profile. 

Open with DEXTER  
In the text 
Fig. 14
The effective magnetic field response over the absolute noise level (logarithmic scale) for the COG method (open circles) and the OMP method (solid circles). The response curve shows how much of the noise is erroneously attributed to a magnetic field. On the left side, the response is shown for a noise range between 10^{5} and 10^{3}. The noise response of the COG method rapidly increases with the noise level from 2 to 100 G, whereas the OMP method remains largely unaffected by the noise with a flat response not larger than 0.7 G. On the right side, a closeup of the magnetic response curve for the region between 10^{5} to 10^{4} is shown. Even in such a low noise regime, the COG method falsely interprets a significant portion of the noise as magnetic field. The OMP method stays essentially flat at a level close to zero over the entire region. The gray shaded area marks the region of incorrect response for the COG method, i.e. , and defines the lower limit for magnetic field detections 

Open with DEXTER  
In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.