Issue 
A&A
Volume 530, June 2011



Article Number  A148  
Number of page(s)  23  
Section  The Sun  
DOI  https://doi.org/10.1051/00046361/201016426  
Published online  27 May 2011 
Validated helioseismic inversions for 3D vector flows^{⋆}
^{1}
MaxPlanckInstitut für Sonnensystemforschung, MaxPlanckStraße 2, 37191 KatlenburgLindau, Germany
email: svanda@mps.mpg.de; gizon@mps.mpg.de; hanasoge@mps.mpg.de
^{2}
Institut für Astrophysik, GeorgAugustUniversität Göttingen, 37077 Göttingen, Germany
^{3}
Department of Geosciences, Princeton University, Princeton, NJ 08544, USA
^{4}
Keldysh Institute of Applied Mathematics, Russian Academy of Sciences, Miusskaya Sq. 4, 125047 Moscow, Russia
Received: 31 December 2010
Accepted: 20 April 2011
Context. According to timedistance helioseismology, information about internal fluid motions is encoded in the travel times of solar waves. The inverse problem consists of inferring threedimensional vector flows from a set of traveltime measurements. While only few tests of the inversions have been done, it is known that the retrieval of the smallamplitude vertical flow velocities is problematic. A thorough study of biases and noise has not been carried out in realistic conditions.
Aims. Here we investigate the potential of timedistance helioseismology to infer threedimensional convective velocities in the nearsurface layers of the Sun. We developed a new subtractive optimally localised averaging (SOLA) code suitable for pipeline pseudoautomatic processing. Compared to its predecessor, the code was improved by accounting for additional constraints in order to get the right answer within a given noise level. The main aim of this study is to validate results obtained by our inversion code.
Methods. We simulate traveltime maps using a snapshot from a numerical simulation of solar convective flows, realistic Born traveltime sensitivity kernels, and a realistic model of traveltime noise. These synthetic travel times are inverted for flows and the results compared with the known input flow field. Additional constraints are implemented in the inversion: crosstalk minimization between flow components and spatial localization of inversion coefficients.
Results. Using modes f, p_{1} through p_{4}, we show that horizontal convective flow velocities can be inferred without bias, at a signaltonoise ratio greater than one in the top 3.5 Mm, provided that observations span at least four days. The vertical component of velocity (v_{z}), if it were to be weak, is more difficult to infer and is seriously affected by crosstalk from horizontal velocity components. We emphasise that this crosstalk must be explicitly minimised in order to retrieve v_{z} in the top 1 Mm. We also show that statistical averaging over many different areas of the Sun allows for reliably measuring of average properties of all three flow components in the top 5.5 Mm of the convection zone.
Key words: Sun: helioseismology / methods: data analysis / Sun: oscillations
Figures 16–28 are available in electronic form at http://www.aanda.org
© ESO, 2011
1. Introduction
The subsurface of the Sun is optically thick, preventing us from directly observing the interior layers. Understanding the properties of the plasma in these regions has consequences for the theories of convection, stability of sunspots, the dynamics of stratified convection, and others. Most current knowledge about convection comes primarily from computational work (e.g., Vögler et al. 2005; Benson et al. 2006; Rempel et al. 2009). Helioseismic inversions of the subsurface flows will play an important role in constraining these theories.
A powerful way of imaging the solar interior is via inferences gathered from studying the statistics of the acoustic and surface gravity waves at the surface. Solar pressure and surface gravity modes are generated randomly by the vigorous turbulence in the upper convection zone. These oscillations are observed in the solar photosphere by measuring Doppler shifts of photospheric absorption lines. Forward modelling allows us to relate anomalies (like flows, thermal hot/cold spots etc.) to changes in helioseismic observables.
The aim of helioseismic inversions is to reveal the structure of the subsurface flows (rotation, meridional circulation, convection), magnetic fields, and to measure deviations in the plasma state parameters (temperature, density, pressure) from a quiet Sun average. In this paper, we focus on travel times (Duvall et al. 1993), i.e., quantities that emerge from fits to cross correlations of observed signals. Timedistance helioseismology is used to measure and interpret changes in travel times of seismic waves caused by inhomogeneities in the structure of the Sun (see review by Gizon et al. 2010). In recent years, timedistance helioseismology has been used to invert for nearsurface flows (e.g., Gizon et al. 2000; Duvall & Gizon 2000; Zhao & Kosovichev 2004; Jackiewicz et al. 2008), for flows beneath sunspots (e.g., Duvall et al. 1996; Zhao et al. 2001; Couvidat et al. 2006; Cameron et al. 2008; Gizon et al. 2009; Moradi et al. 2010) and flows in their vicinity (Gizon et al. 2000), study the rotational gradient at the base of the convection zone (e.g., Hanasoge & Duvall 2009), etc.
Helioseismic inversions are performed using two principal methods: the regularised least squares (RLS) and optimally localised averaging (OLA). The RLS method (in timedistance helioseismology used for the first time by Kosovichev 1996) seeks to find the models of the solar interior, which provide the best leastsquares fit to the measured traveltime maps, while regularising the solution (e.g., by requiring the smooth solution). The OLA method was developed for geoseismology (Backus & Gilbert 1968, 1970). A form suitable for use in helioseismology was devised by Pijpers & Thompson (1992), who formulated the subtractiveOLA method. SOLA is based on explicitly constructed spatially confined averaging kernels by taking linear combination of sensitivity kernels, while simultaneously keeping the error magnification small. The resulting coefficients are then used to linearly combine the traveltime maps and obtain the estimate for structure and magnitude of solar plasma perturbations. A SOLAtype inversion is the principal method discussed in the current paper. The SOLA has been used in timedistance local helioseismology in the past by Jackiewicz et al. (2007, 2008) who demonstrated the ability of SOLA inversions to reveal the structure of the 3D internal flows. An efficient approach to solve fully consistent SOLA inversions was introduced by Jackiewicz et al. (2011). In this paper we focus on inversions for threedimensional vector flows on supergranular scales in the nearsurface layers of the solar convection zone.
1.1. Validating helioseismic inversions
Two approaches have been used so far to validate timedistance inversions. The first approach consists of generating synthetic traveltime maps by convolving a (known) arbitrary frozen flow field with traveltime sensitivity kernels and then testing the inversion method using these synthetic travel times (e.g., Kosovichev & Duvall 1997). The second approach (Zhao et al. 2007, 2010) is to use evolving realistic numerical simulations of threedimensional radiative convection, where the helioseismic waves are naturally excited by the convection.
The first approach is convenient, but may not represent a realistic situation, in particular because the prescribed subsurface structures are often too idealised. The second approach is preferable, but is limited by computing resources: only simulations of the very nearsurface regions of the Sun are available today. Both approaches have shown that inversions are generally able to retrieve the horizontal components of velocity at supergranulation scales in the quiet Sun. However, some problems have been reported. For example, vertical velocities have been measured with the opposite sign near the surface (Zhao et al. 2007) and the structure and sign of the flows around sunspots vary with the inversion method used (see, e.g., Gizon et al. 2009; Moradi et al. 2010).
The aim of this paper is to validate a particular implementation of SOLA inversions for timedistance helioseismology using a mixed approach: we take a snapshot from a largebox realistic simulation of solar convection (Ustyugov 2008) to generate realistic traveltime maps by convolution with Born traveltime sensitivity kernels. A realistic noise component is added to the travel times. The traveltime maps are then inverted using a multichannel SOLA inversion and compared with the known flows. This approach allows us to investigate various types of biases in the results of the inversion and to develop a robust procedure to minimise them. The most serious bias comes from the natural correlations among the components of flows induced by mass conservation that translate into a crosstalk between the components of the inverted flow (as we shall see, this is important when retrieving the vertical flow). Using realistic numerical simulations of convective flow velocities as input is very useful to set the acceptable level of random noise of the inverted flow velocities. The use of travel times with realistic noise properties is important in order to derive realistic estimates of the noise of the inverted flow velocities.
2. Synthetic travel times
To construct synthetic travel times, we use a realistic hydrodynamic simulation (Ustyugov 2008) of the solar convection. The computational domain is a box 20 Mm in depth and 60 Mm in each horizontal direction. The simulation provides us with a reasonably realistic description of flows in the upper convection zone v = (v_{x},v_{y},v_{z}). Throughout this paper the spatial coordinates are defined as (1)where r is the horizontal position vector and z is the height.
We choose a snapshot (t = 500 min, see Fig. 16 in the online supplement) of the above mentioned simulation for our inversion tests. Following Duvall et al. (1997), we consider different types of traveltime measurements between a surface point at position r and a concentric annulus or quadrants in order to measure travel times sensitive to flows in inwardoutward (denoted by “oi”), westeast (“we”), and northsouth (“ns”) directions. Traveltime maps are denoted by τ^{a}(r), where the superscript a is an integer that uniquely refers to a particular combination of choices in the data analysis: a type of geometry (oi, we, or ns), annulus radius (from 7.3 Mm to 29.2 Mm every 1.46 Mm), and a wave filter (here, ridge filters for one of f, p_{1}, p_{2}, p_{3}, or p_{4} modes). Thus, index a refers to one of M = 3 × 16 × 5 = 240 possibilities. Additional information describing the measurement procedure in timedistance helioseismology is given by Gizon & Birch (2005).
Traveltime maps (travel times as functions of position vector r) are generated by convolving the convection snapshot with sensitivity kernels according to (2)where is a vector traveltime sensitivity kernel (see, e.g., Fig. 17 in the online supplement), v is the velocity vector of convecting flows, and the volume integral is taken over the Sun. The noise component is denoted by n^{a}. The original horizontal size of the simulation box of 60 Mm was too small for our purpose. Since the simulated velocities are periodic in the horizontal directions, we copied the simulation box 10 × 10 times.
Kernels are computed using the singlescattering Born approximation (Birch & Gizon 2007) and depend on eigenmodes of a background 1D standard solar model (ChristensenDalsgaard et al. 1996). The kernels are invariant under horizontal translations; horizontal averages of these kernels for various oscillation modes are displayed in Fig. 18 in the online supplement to display their sensitivity in depth. All kernels used in this study have sensitivities only in uppermost 10 Mm of the convection zone.
Solar waves are excited by the action of turbulent convection in the Sun and therefore travel times are inherently noisy. Following Gizon & Birch (2004), we use a realistic noise covariance matrix (see Fig. 19 in the online supplement), (3)to generate realizations of the noise.
Each traveltime map has 400 × 400 pixels with spatial sampling of 1.46 Mm, corresponding to the pixel size of a MichelsonDopplerImager fulldisc image (MDI, Scherrer et al. 1995). In total, we generate 240 different traveltime maps (one for each index a, see, e.g., Fig. 1). These traveltime maps have spatial power spectra that are similar to observed traveltime power spectra based on MDI data, as demonstrated in Fig. 2. However, our synthetic travel times have weaker power at low wavenumbers due to the complete lack of signal from the simulated convective velocities below k R_{⊙} = 60 (only simulated noise is present).
Fig. 1 Examples of synthetic traveltime maps used in this paper (top row) and observed travel times from SOHO/MDI (bottom row). Temporal length of observation is T = 6 h. 
Fig. 2 Comparison of azimuthally averaged spatial power spectra of synthetic (solid) and observed (dashed) travel times displayed in Fig. 1 as a function of k R_{⊙}, where k is the horizontal wavenumber and R_{⊙} is the radius of the Sun. 
We check to ensure that pointtoannulus sensitivity kernels and noise covariance matrices obey expected symmetries. Furthermore, the sensitivity kernels and covariance matrices are forced to decay smoothly to zero towards the edge of the computation box by multiplying with a smooth spatial function having zeros far from the region of interest. In addition, we make sure that the horizontal integral of each sensitivity kernel for v_{z} is zero at all depths. This is a consequence of symmetries associated with pointtopoint kernels for v_{z}, which are insensitive to the mean value of vertical velocity (Birch & Gizon 2007).
Notations.
3. Improved SOLA inversions
3.1. Realspace formulation
The inversion method we use here is an improved version of the subtractive optimally localised averages (SOLA) method described by Jackiewicz et al. (2007, 2008, 2011).
The SOLA algorithm describes how to optimally combine a given set of travel time measurements to infer the underlying properties of medium. In this paper, we wish to retrieve , i.e., an estimate of the αcomponent of the flow velocity in the neighbourhood of position x_{0} = (r_{0},z_{0}) in the solar interior. In practice, we search for a linear combination of the traveltime measurements, (4)where are inversion weights to be determined. Combining Eqs. (2) and (4), we have (5)where is a vector averaging kernel with components (6)The expression d^{2}r dz means d^{3}x.
The expectation value of inverted velocity may be expanded as (7)where (8)and . The form of Eq. (7) will be useful later in this paper to optimize the inversion. Notice that the second term on the right side of (7) represents leakage of other flow components into the inverted one. This crosstalk may be a significant source of bias and needs to be studied.
We search for inversion weights w^{α} so that the vector averaging kernel resembles a usersupplied target . Because we wish to invert for the αcomponent of the flow, we choose the target with components: (9)where is a function that peaks around x_{0} = (r_{0},z_{0}), diminishes rapidly away from that point, and has unit integral. δ_{αβ} is a Kronecker δ. Throughout this paper, we choose a simple Gaussian (10)where s_{h} and s_{z} are full widths at half maximum (FWHM) of the Gaussian in horizontal and vertical directions respectively.
A successful inversion will return a value of close to the target velocity (11)and crosstalk velocities close to zero for β ≠ α.
A list of the notations used for the various flow velocities referred to in the inversion procedure is given in Table 1.
The problem to be solved to obtain the w^{α} belongs to the class of constrained regularised optimisations. The terms to be regularised are the following:

The misfit, i.e., how far the averaging kernel is from the desired target function, (12)where is an implicit function of the weights w^{α}.
The crosstalk quantifying the leakage of the signal from other flow components into the inverted one, (13)

The variance of , which corresponds to the rootmeansquare of noiserelated fluctuations in inferred velocities. We refer to this quantity as predicted error hereafter, (14)

The ad hoc term quantifying the spread of the weights in space around r_{0}, (15)The regularisation based on S ensures that weights decrease towards the edge of the horizontal domain in order to enforce the spatial locality of the inversion and to prevent the weights from oscillating in the spatial domain.
The inversion is subject to constraints (17)in order to scale the amplitude of the inverted flow appropriately.
By taking derivative of (16) with constraint (17) added with respect to the weights, the problem can be cast into a linear inverse problem, as explained by Pijpers & Thompson (1992) or Jackiewicz et al. (2011). In real space, the matrix to be inverted has (N^{2}M + P)^{2} ≃ (10^{7})^{2} elements, where N = 200 is the number of grid points in one horizontal direction, P = 3 is the number of physical unknowns (three velocity components), and M = 240 as already defined earlier.
3.2. Fourierspace (multichannel) formulation
The full problem written in real space is intractable to be solved using nowadays computers, because the matrix to be inverted is too large. Jackiewicz et al. (2011) found a solution to this problem by transforming to spatial Fourier space, where the inverse problem decouples as a consequence of the horizontal translation invariance of the sensitivity kernels. Thus, instead of a big linear inverse problem (inverting matrix having (10^{7})^{2} elements), we solve 40 000 small linear inverse problems (inverting matrices having 240^{2} elements) in wavevector space. This approach is called a multichannel inversion (Jensen et al. 1998). Here we summarise results of Jackiewicz et al. (2011) for the sake of completeness.
Following Jackiewicz et al. (2011), we use the following definition of the Fourier transform, such that any function f(r) and its 2D spatial Fourier transform are related according to where h_{x} = 1.46 Mm and h_{k} = 0.022 rad/Mm are the grid spacings in the real and Fourier domains respectively.
For each nonvanishing wavevector k = (k_{x},k_{y}), the vector of weights is the solution to the matrix equation, (20)Each matrix A(k) has M × M elements A_{ab} given by (21)The vector has M elements given by (22)In addition, the matrix equation for the case k = 0 is (23)where C is an M × P matrix with elements (24)L is a 1 × P vector of Lagrange multipliers, and the vector U^{α} is a 1 × P unit vector with components .
Solutions to the N − 1 Eq. (20) and the Eq. (23) give Fourier components of the weights, . By taking an inverse Fourier transform we obtain . Equation (4) then gives an estimate of the αcomponent of the flow, .
The difference between the above equations and those in Jackiewicz et al. (2011) arises from additional constraint terms (13) and (15).
3.3. Picking tradeoff parameters
Tradeoff parameters μ, ν, and ϵ control the balance between various terms in the costfunction (16).
In practice, during the inversion for each flow component, we compute a grid of solutions by varying all three tradeoff parameters. For each solution on this grid, we compute the misfit (12), predicted error of results (14), amount of crosstalk (13), and the spatial power of weights (15). All four quantities can be computed in spatial Fourier domain.
Standard optimisation methods to achieve optimality, such as Lcurve analysis (see, e.g., Hansen 1998; Jackiewicz et al. 2008), are not particularly useful to our problem. The elbow of the Lcurve, which is considered an optimal point in parameter space, was located where the predicted noise level was very large. We therefore developed our own strategy of selecting values for the tradeoff parameters.
We start with parameter μ, which controls the tradeoff between the misfit and random error. We select μ such that the noise level of the inverted quantity is less than some target value, chosen by the user, depending on the problem at hand. For example, one may choose a target noise level of 20 m/s to invert for horizontal flows in supergranules. A much lower noise level will be required to invert for the vertical component of the flow (Duvall & Birch 2010). The example noise levels σ_{α} are given in Table 2, where we give the predicted inversion error for and and the rootmeansquare of the flow averaged with the target function ( and ) and averaged with the resulting averaging kernel ( and ). It may be that this selection of μ leads to an averaging kernel which does not resemble the desired target function at all, in which case one would have to allow for a coarser spatial resolution.
Predicted inversion errors and the expected magnitude of the velocities at the three different depths discussed throughout this paper.
For the given μ, we then choose ν so that the degree of crosstalk is less than 10^{5} Mm^{3}. We observe that this constraint places the following upper bound on the magnitude of the crosstalk terms: (25)The parameter ϵ controls the degree of spatial confinement of the inversion weights and has an impact on the misfit and error as well. If ϵ is too large, the weights will be highly localised in space around the central point, but a large misfit will result.
For the particular cases that we study here, we find that localisation of weights is accomplished when S < 2 × 10^{3} km^{2} s^{2}. The impact of this regularisation term is illustrated in Fig. 3. For the case shown, weights are confined to within a disc of radius ~50 Mm.
Fig. 3 Example inversion weights for 1 Mm depth and effect of the minimisation of the spread S. A strong regularisation ϵ confines the inversion weights in the spatial domain (bottom row, ϵ = 10^{0.5}) compared to the case ϵ = 0 (top row). The cut at y = 0 (horizontal black line) is displayed on the righthandside panel, where the values of were scaled by . 
3.4. The code
The code implementing the above described procedure is written in Matlab. Matlab provides a compromise between computational efficiency and the availability of higherlevel software constructs that make the code lucid, modular, and easy to modify. Matlabbased code may be compiled into binaryexecutable form suitable for pipeline pseudoautomatic processing. The input is a text file, in which the user specifies all relevant parameters including the kernels to be used in the inversion and the set of tradeoff parameters to be investigated. This allows for the development of user interfaces such as the webbased or graphical interface, which can serve as the user frontend of the inversion code.
The code is parallel and scales linearly with number of parallel jobs. The execution is fast, e.g., a set of inversions for one flow component involving 240 kernels, 200 × 200 spatial points, and a grid of 200 tradeoff parameter values takes around 6 h using 48 Opteron2.3 GHz CPUs. When many CPUs are involved, the extensive inputoutput load becomes a bottleneck and affects the total processing time.
4. Inversion for horizontal flow
We applied the abovedescribed method to known synthetic travel times in order to validate the setup and the performance of the inversion procedure. We focus first on the horizontal flow components, i.e., α = x or α = y. We compare two principal quantities: (1) the flow map that is obtained by convolving the known velocity field with the target function (), i.e., the bestcase inversion scenario, and (2) the flow map actually resulting from the inversion (). In this manner, we may investigate in detail different sources of bias in the results.
The simulation convolved with the target function gives us , from which we estimate the expected magnitude of the flow that we want to invert for (see Table 2). This places limits on the required noise level of the inversion so that the results have signaltonoise ratios larger than 1. Requirements on the targeted error level fix the tradeoff parameter μ.
The choice of the target function depends on the discretion and needs of the user. Here we focus on layers in the top few Mm, in particular on depths − z_{0} = 1, 3.5, and 5.5 Mm. This setup was selected because similar flow inversions were also performed and discussed by Jackiewicz et al. (2008), which makes it possible to compare the results of both methods. As discussed in the preceding sections, the outcome of the inversion is a set of weights , which are used to combine the traveltime maps in order to obtain estimates of velocity and . Example weights for are displayed in Fig. 3.
All components of the averaging kernels are shown in Figs. 20–22. We also show comparisons between inversions when the crosstalk is minimised and not minimised. The minimisation of the crosstalk with ν = 100 is very efficient. However, for depths of 3.5 and 5.5 Mm the minimisation of the crosstalk introduces some small artefacts in the inversion averaging kernel , which are a small price to pay.
Figure 4 shows vertical cuts through the target functions and the averaging kernels when the crosstalk is minimised. Except for the target depth of 5.5 Mm, the depth dependences of the averaging kernels resemble that of the target functions. The nearsurface inversion at 1 Mm is dominated by the fmode. The inversion at 3.5 Mm has equal contributions from f and p_{1}, with opposite signs, as shown in Fig. 23 in the online supplement, where we plot the contribution of individual modes to the averaging kernels as a function of depth.
Fig. 4 The cut through the x = y = 0 point of the averaging kernel (solid) and the respective target function (dashed) for v_{x} inversion with minimised crosstalk, using travel times averaged over 4 days, at three discussed depths (1, 3.5, and 5.5 Mm). 
Fig. 5 Comparison of inverted v_{x} with input data. Toprow panels show the input flow field convolved with the target function and the averaging kernel respectively. Bottom rows show inversion results in cases when crosstalk is ignored (left) and minimised (right). Random errors of the inversion are given in Table 2. 
Fig. 6 Inversion biases for v_{x} at 1 Mm depth. Left: noiseless with crosstalk minimised versus the ideal . The departure from slope unity (grey dashed line) is due to an imperfect match between and . Right: noisy versus noiseless . The results are plotted in two cases: when the crosstalk is ignored (grey × ) and minimised (black + ). The linear fit to the black crosses coincides with the dashed line of slope unity. The black dotted lines represent the predicted error of 14 m s^{1}, which is consistent with the observed scatter of the black crosses. 
Fig. 7 The azimuthallyaveraged power spectra of various components of the v_{x} inversions using travel times averaged over 4 days. For reference, we plot the power spectrum of (thick grey solid line). We plot also the power spectrum of (solid line) and the power spectrum of the noise (i.e., the power spectrum of ; dashed line) for the inversion where the crosstalk is minimised (black) and ignored (grey). 
The validation of the v_{x} inversion is demonstrated in Fig. 5. Here we plot the desired at three different depths and the inverted without noise contributions. These two are very close for the depths 1 and 3.5 Mm. The differences between and at depth 5.5 Mm are caused by an imperfect averaging kernel. When the random noise is added to the solution (bottom row of Fig. 5), we see that the inversions for v_{x} at 1 Mm and 3.5 Mm are still very good, while the inverted v_{x} at depth 5.5 Mm is dominated by noise. Minimising crosstalk does not improve the quality of the solution in this case.
In order to quantify the inversion biases, we compare directly the expected and inverted values for a set of spatial locations. We select points separated by 7.5 Mm in each horizontal coordinate. This sampling interval is equal to the half of the horizontal FWHM of the target function (thus the points are somewhat independent). In addition, to avoid possible edge effects, we cut the outer part of the horizontal plane so that only the central 200 × 200 pixel patch is kept. The scatter plots comparing various inversion components are displayed in Fig. 6 for the depth of 1 Mm, the plots for the depth of 3.5 Mm are qualitatively similar. We estimate that for the depth of 1 Mm, the imperfect averaging kernel leads to an average underestimation of the horizontal flow components by some 20% (Fig. 6 left). Furthermore we note that in the case of the inversion for the horizontal flow components, the bias caused by the crosstalk is not important and the random noise level corresponds to the predicted value (Fig. 6 right). The results for the depth 3.5 Mm are similar, while the results for the depth 5.5 Mm are dominated by the random noise.
Correlation analysis of various inversion components using synthetic travel times averaged over 4 days.
The noise domination of v_{x} inversion at depth of 5.5 Mm is evident from the plot of the azimuthally averaged power spectra of the random noise and of the signal as a function of k R_{⊙} (Fig. 7 right). For comparison, similar plots for the depths 1 Mm and 3.5 Mm are given also in Fig. 7, where the power spectrum of the signal is well above the random noise near supergranular spatial scales. The decrease in power of the signal at low kR_{⊙} in Fig. 7 is because the convection simulation does not contain these scales.
We showed that it is possible to retrieve v_{x} and v_{y} in the top 3.5 Mm without noticeable bias and within the predicted noise level of ~25 m/s (for observing time T = 4 days) with 240 different traveltime measurements for ridges f to p_{4}. In Table 3, we summarise our findings. We state the statistical quantities (correlation coefficient and the slope of the linear fit) comparing and the corresponding . We also show the signaltonoise ratio of the inverted .
The crosstalk is unimportant in inversions for the horizontal flow components. The crosstalk could come only from v_{z}, which is weak. We estimate that the influence of the crosstalk is less than 5%. The inversion at 5.5 Mm is already dominated by noise. If we drop the requirement on horizontal resolution, inversions at a depth of 5.5 Mm are also possible (e.g., with a FWHM of s_{h} = 25 Mm).
It is possible to perform inversions with less temporal averaging, but to obtain a reasonable signaltonoise ratio, we would have to relax the demand on the match between the averaging kernel and the target function. Typically, this leads to sidelobes in the averaging kernel (especially in the zdirection), which may make interpretation of the results more difficult. Another possibility would be to include a larger number of independent traveltime measurements (and therefore more sensitivity kernels). These issues are being worked on.
5. Inversion for vertical flow
5.1. Specificity of inversions for vertical flow
Inversions for v_{z} require a different methodology because sensitivity kernels have zero horizontal integrals at each depth (and therefore zero total integral) implying (see Eq.(6)) (26)Consequently, with measurements discussed here, it is impossible to retrieve horizontal average ⟨ v_{z} ⟩ of vertical flow and we may only invert for fluctuations v_{z} − ⟨ v_{z} ⟩ . Equation (26) implies
(27)and thus Eq. (17) cannot be written for α = z. This implies in turn that the matrix Eq. (23) for k = 0 may not be written either.
5.2. Ignoring k=0
A solution to this issue is simply to replace Eq. (23) by (28)Since the averaging kernel must have zero integral, the averaging kernel will be offset by a small negative constant away from the central peak. Luckily, regularization by the term (15) will ensure that the averaging kernel will drop to zero at the edge of the inversion box. However, this solution is not quite satisfactory, because we do not have very much control over this extended negative surrounding sidelobe.
5.3. Target function with zero mean
It is more elegant to select a target function with vanishing horizontal integral at each depth. If this target function is to peak around x_{0} = (r_{0},z_{0}), then this peak must be compensated by negative sidelobes in horizontal directions. These sidelobes have to be constructed in such a way that they do not lead to significant biases. We suggest to choose (29)where (30)The constant (31)ensures that the horizontal integral of is zero. The horizontal part of the target function is constructed from a Gaussian peaked at r_{0} = 0, from which a wide surrounding Gaussian annulus is subtracted. The free parameter n balances the width of a sidelobe compared to the width of the central peak.
The inversion is performed for each k ≠ 0 with this new target function, together with .
5.4. Both target functions provide similar answers
We compared the results obtained using the two target functions proposed above and found two conclusions.
The two solutions described above provide results which are very close (see, e.g., Fig. 24). However, the solution with the target function having zero horizontal average is elegant and provides more control over the solution averaging kernel.
By comparing the results obtained with various sizes of the negative annulus surrounding the central peak we found that the bias in the inverted flow caused by this negative sidelobe is negligible if n > 3 in Eq. (30).
The results presented in the following sections are obtained using the first formalism.
Fig. 8 Horizontal spatial crosscorrelation between v_{x} and v_{z} from the simulation in the top 1 Mm. These flow correlations together with a nonvanishing crosstalk averaging kernel may lead to a systematic bias of the inferred v_{z}. The circle of radius of 15 Mm shows the width of the desired target function. 
Fig. 9 The cut through the x = y = 0 point of the averaging kernel (solid) and the respective target function (dashed) for v_{z} inversion using travel times averaged over 4 days. Compare with Fig. 4. 
5.5. Validation of vertical flow inversion
As is evident already, inversions for vertical flow are not as seamless as in the case of horizontal components. This is mostly because the vertical flow is much weaker on supergranular scales than the horizontal flow. As shown in Table 2, the expected rms of the vertical flow in the top layers is of the order of 5 m s^{1}. Therefore, the predicted noise of the inverted v_{z} has to be set to a much smaller value than for the horizontal velocities. As a result, the match between the desired target function and the resulting averaging kernel will have to be poorer.
Furthermore, it is absolutely crucial to minimise the crosstalk to avoid the leakage of the large horizontal velocities into the small inferred vertical velocity. Even an apparently negligible crosstalk averaging kernel could in the end cause a significant bias in the results. Minimising the crosstalk is especially important because of the natural correlations between the vertical and horizontal flow components in the massconserving flow of the supergranules (see Fig. 8). In the upper layers, horizontal outflows are associated with upflows. The structure of the crosstalk averaging kernel as shown in the top half of Fig. 10 in the case when crosstalk is not minimised then implies a negative bias in v_{z} due to and .
Fig. 10 All components of the averaging kernel for v_{z} inversion at 1 Mm depth with a FWHM of s_{z} = 1.1 Mm and s_{h} = 15 Mm. Bottom row: with crosstalk minimised, top row: crosstalk is ignored. The crosstalk is presented in the form of and averaging kernel components. Random error of the results is 3 m s^{1} when assuming data averaged over 4 days. Overplotted contours, which are also marked on the colour bar for reference, denote the following: halfmaximum of the kernel (white), halfmaximum of the target function (black), and by grey lines ± 5% of the maximum value of the kernel (solid and dotted, respectively). 
Fig. 12 All components of the v_{z} inversion at 1 Mm depth. Top row with the crosstalk ignored, bottom row with the crosstalk minimised. We demonstrate that if crosstalk is not addressed, horizontal components will leak into the inverted v_{z} and cause a bias. 
Figure 9 shows vertical cuts through the target functions and the averaging kernels when the crosstalk is minimised. The averaging kernels for the inversions at depth 3.5 and 5.5 Mm have sidelobes towards the surface. All components of the averaging kernel for v_{z} are displayed in Fig. 10 (and Figs. 25, 26), again comparing the cases when the crosstalk is and is not minimised. The action of the crosstalk minimisation term is very efficient.
Fig. 13 Inversion biases for v_{z} at 1 Mm depth, similar to Fig. 6. While the black line (with minimised crosstalk) almost coincides with the dashed line of slope unity, the grey one (crosstalk ignored) indicates a horrible bias. The predicted error is 3 m s^{1}. The inversion is not sensitive to the horizontally averaged horizontal flow ⟨v_{z}⟩ = 11 m s^{1}, therefore it was subtracted from in this comparison. 
Fig. 14 The azimuthallyaveraged power spectra of various components of the v_{z} inversions using travel times averaged over 4 days. Compare with Fig. 7 for v_{x} inversions. Notice that in the case of containing the crosstalk at 1 Mm depth, the excess in power around kR_{⊙} ~ 150 is not caused by the random noise, but by the bias coming from the crosstalk (not plotted separately). In the other two cases the random noise is a main cause why these inversions are not possible. 
The validation of the v_{z} inversion is demonstrated in Fig. 11. Here we plot at three different depths and the inverted without noise contributions. These two are very similiar at the depth of 1 Mm. Small differences between and caused by an imperfect averaging kernel are visible at the depth of 3.5 Mm, become more significant at depth 5.5 Mm. When random noise is added to the solution (bottom row of Fig. 11), we see that the inversion for vertical flow is possible at 1 Mm depth only if crosstalk is minimised. Inversions at depths 3.5 and 5.5 Mm are buried in random noise, where minimising crosstalk does not help. The magnitude of crosstalk at a depth 1 Mm is much larger than that of the vertical flow; further, crosstalk is highly anticorrelated with the vertical flow, which makes the correlation coefficient between and close to − 1. Note that a similarly high anticorrelation was measured by Zhao et al. (2007).
All components and the random noise component in the inversion for v_{z} at 1 Mm depth are displayed in Fig. 12. We see in the top row that the leakage of the horizontal components ( and ) covers up completely the weak signal of the vertical flow when the crosstalk is not minimised.
Biases in the v_{z} inversion at a depth of 1 Mm may be quantified by directly comparing expected and inverted values for a set of spatial locations in the horizontal plane (Fig. 13), as in the case of the v_{x} inversion (Sect. 4). The imperfect averaging kernel in this case does not cause any significant bias (Fig. 13 left). The v_{z} inversion is not sensitive to the horizontally averaged vertical flow ⟨ v_{z} ⟩ , which is therefore subtracted from in the corresponding plot. The effect of minimising crosstalk is shown in Fig. 13 left. The inverted is anticorrelated with expected from the inversion due to the leakage of the horizontal massconserving flow components into the vertical one.
Meaningful inversions for the vertical flow are possible in shallow nearsubsurface layers of the convection zone only when the crosstalk between the vertical and the horizontal components is minimised. Vertical flow inversions on supergranular horizontal scales at depths greater than ~1 Mm performed using f to p_{4} modes and travel time maps averaged over 4 days are dominated by random noise and the signal of the vertical flow may not therefore be inferred at all (see summary in Table 3). This fact is also demonstrated in Fig. 14, where we plot power spectra of individual inversion components as a function of the spatial scale. At depths larger than 1 Mm, the power of the signal is much less than the power of the random noise regardless of the spatial scale.
Fig. 15 Comparison at three depths of statistical averages of velocityvector over many flow realisations. The horizontal flow is displayed by arrows while the vertical flow is colourcoded. The reference arrow indicates 100 m s^{1}, the random errors of the inversions are given in Table 2. Compare to Figs. 5 and 11. 
6. Beating the noise: statistical averaging
As demonstrated in previous sections, it is very difficult to have a meaningful inversion for solar flows on supergranular scales averaged over a short time, even in the shallow nearsubsurface layers of the solar convection zone. The signal is overwhelmed by random noise. It is not feasible to average over longer times, because the expected lifetime of convection on supergranular scales is on the order of a day. By averaging over longer time, a significant portion of the scientifically useful information would be lost.
It seems feasible to solve the issue by averaging over many realisations of similar flow structures each averaged over short time (therefore noisy). This concept was already used by other authors (e.g., Duvall et al. 2006). As an example we note the possibility to measure the flows in many individual supergranules and to average over this sample in order to obtain the typical flow structure in an “average supergranule”, similar to Duvall & Birch (2010). The predicted error in the results after averaging scales as , where is the size of the ensemble of independent flow realisations. This scaling law allows us to relax the noise constraint on the inversion and therefore regularise more strongly about the misfit and crosstalk terms. See Table 2 for required random error of inversions results.
For instance, let us assume that we average over 10^{4} independent inverted flows, each obtained by considering traveltime maps averaged over 6 h. Individual flow maps are noisy and therefore contain little useful information about underlying flows. This selection allows us to relax the constraint on the estimated error level by a factor of a 100 and therefore obtain a much better fit to the target function (see Fig. 27). With this setup, the validation displayed in Fig. 15 shows almost perfect correspondence between expected v^{tgt} and inverted flows in the top 5.5 Mm of the convection zone. Inversion for v_{z} at 5.5 Mm with a FWHM of s_{h} = 15 Mm is an exception, because the selected target function cannot be matched by the averaging kernel. Different inversions for v_{z} at this depth are possible with greater averaging, e.g., with s_{h} = 25 Mm. It is still crucial to minimise crosstalk in the case of v_{z} inversions to retrieve the correct answer. Power spectra of inversion components at all discussed depths are displayed in Fig. 28.
7. Conclusions
We improved and validated an inversion algorithm based on a SOLA inversion approach. The formalism, algorithm and the code is universal and can, in principle, be used to invert for any quantity describing inhomogeneities in the solar plasma, provided that the corresponding sensitivity kernels are available. The code is also ready to be used for application to real measurements. It will become part of the helioseismic pipeline running at German Data Center for SDO. We plan to use it to analyse all available SOHO/MDI and SDO/HMI data in order to routinely provide a tomographic image of the structure of the solar upper convection zone.
The code is absolutely scalable allowing to include more independent measurements (and therefore more sensitivity kernels) in order to further refine the precision of the results. Thanks to the decoupling, the problem will still be solvable even using nowadays computers.
We improved the inversions by introducing additional terms that allow to control and minimise some sources of bias in the results. Most importantly, we minimise the crosstalk between individual flow components, which is crucial especially for v_{z} inversions. In principle, the formalism allows to minimise crosstalk between any selected quantities during the general inversion. The validation performed here silently assumes that the sensitivity kernels and noise covariance matrices are perfect. For the validation of the method and the code it is important if the sensitivity kernels and noise covariance matrices are solarlike. This allows us to study different sources of biases we may expect in the Sun. It is also very important that the sensitivity kernels used in the inversion contain all the details of the traveltime measurements, including the instrumental function affecting the solar oscillation power spectrum.
We found that by considering f to p_{4} frequencyaveraged modes and supergranular spatial scales it is possible to perform reliable and trustworthy flow inversions of the traveltime maps averaged over a few days in the top 3.5 Mm layer of the convection zone in the case of horizontal v_{x} and v_{y} components, and in the top 1 Mm in the case of vertical v_{z} component. Based on our experiment, we expect that using traveltime maps averaged over 4 days, it is possible to measure 3D velocities as weak as 10 m s^{1} at the surface and horizontal velocities having amplitude 20 m s^{1} in the top 5 Mm. We estimate, that using the traveltimes averaged over 1 day, it should still be possible to measure all components of the supergranular flow at the surface and its horizontal components in the top 5 Mm. The crosstalk minimisation is crucial in order to measure the correct vertical velocity. Its presence in the results may explain the opposite sign of the vertical flow inversion discovered by Zhao et al. (2007).
By considering many flow realisations and statistical averaging we might go deeper in the convection zone to learn about the horizontal flow components and perform vertical flow inversions for the top few Mm depths.
Online material
Fig. 16 Slices through the x and z components of the simulated flow v at depths 1 Mm, 3.5 Mm, and 5.5 Mm. 
Fig. 17 The pointtoannulus sensitivity kernels for three flow components computed for the fmode, distance 28 Mm and outwardinward geometry. The white circle represents the location of the averaging annulus. 
Fig. 18 Horizontal averages of sensitivity kernels may serve useful when estimating which depths are easier to target. The trend for each mode/ridge was obtained by taking and averaging over all as within the given mode. 
Fig. 19 An example noisecovariance matrix for fmode travel times averaged over 6 h. In this plot, a stands for the combination of the fmode, oi geometry, and annulus radius of 7.3 Mm, b stands for the combination of fmode, we geometry, and annulus radius of 8.8 Mm. 
Fig. 20 All components of the averaging kernel for v_{x} inversion at 1 Mm depth with a FWHM of s_{z} = 1.1 Mm and s_{h} = 15 Mm. Bottom row: with crosstalk minimised, top row: crosstalk is ignored. Random error of the results is 14 m s^{1} when assuming data averaged over 4 days. Overplotted contours, which are also marked on the colour bar for reference, denote the following: halfmaximum of the kernel (white), halfmaximum of the target function (red), and ± 5% of the maximum value of the kernel (blue and green, respectively). 
Fig. 21 All components of the averaging kernel for v_{x} inversion at 3.5 Mm depth with a FWHM of s_{z} = 2.2 Mm and s_{h} = 15 Mm. Random error of the results is 20 m s^{1} when assuming data averaged over 4 days. For details see Fig. 20. 
Fig. 22 All components of the averaging kernel for v_{x} inversion at 5.5 Mm depth with a FWHM of s_{z} = 3.5 Mm and s_{h} = 15 Mm. Random error of the results is 28 m s^{1} when assuming data averaged over 4 days. For details see Fig. 20. 
Fig. 23 The contributions of particular modes to the horizontally averaged averaging kernel for v_{x} inversions using travel times averaged over 4 days for depths 1 and 3.5 Mm. We do not display the inversion for the depth of 5.5 Mm, because it is heavily dominated by noise. 
Fig. 24 To solve the peculiarity of the v_{z} inversion, we introduced two formalisms in Sect. 5.1. Here we plot performance of those. In black, the magnified section a y = 0 and z = z_{0} of different target functions are displayed, one with removed mean (1) and one constructed with negative sidelobes (2). The resulting averaging kernels are also plotted. It is evident that the resulting averaging kernels are qualitatively very similar even when different formalisms were used to compute them. 
Fig. 25 All components of the averaging kernel for v_{z} inversion at 3.5 Mm depth with a FWHM of s_{z} = 2.2 Mm and s_{h} = 15 Mm. Random error of the results is 13 m s^{1} when assuming data averaged over 4 days. For details see Fig. 20. 
Fig. 26 All components of the averaging kernel for v_{z} inversion at 5.5 Mm depth with a FWHM of s_{z} = 3.5 Mm and s_{h} = 15 Mm. Random error of the results is 133 m s^{1} when assuming data averaged over 4 days. For details see Fig. 20. 
Fig. 27 The cut through the x = y = 0 point of the averaging kernel (solid) and the respective target function (dashed) for the v_{x} (left) and v_{z} (right) inversions using averaging over many flow realisations plotted along with the corresponding target functions at three discussed depths (1 Mm in blue, 3.5 Mm in green, and 5.5 Mm in red). Compare to Figs. 4 and 9 where the resemblance of the target functions is worse. The random error of the results is given in Table 2. 
Fig. 28 The azimuthallyaveraged power spectra of the v_{z} inversion components at depths of 1, 3.5 and 5.5 Mm for averaging over many flow representations. For reference, we plot the power spectrum of using the black solid line. Then we plot the power spectrum of (solid line) and power spectrum of the noise (i.e., the power spectrum of ; dashed line) for the inversion where the crosstalk is minimised (blue) and ignored (red). Compare to Fig. 14. 
Acknowledgments
This study was supported by the European Research Council under the European Community’s Seventh Framework Programme (FP7/2007–2013)/ERC grant agreement #210949, “Seismic Imaging of the Solar Interior”, to PI L. Gizon (Milestone #5). Authors would like to thank J. Jackiewicz for providing us with the noise covariance matrices. The flow sensitivity kernels were computed using the code written by A. C. Birch deployed within the HELAS project at http://www.mps.mpg.de/projects/seismo/NA4. M.Š. acknowledges a partial support through the Grant Agency of Academy of Sciences of the Czech Republic under grant IAA30030808.
References
 Backus, G. E., & Gilbert, J. F. 1968, Geophys. J., 16, 169 [NASA ADS] [CrossRef] [Google Scholar]
 Backus, G. E., & Gilbert, J. F. 1970, Roy. Soc. London Philos. Trans. Ser. A, 266, 123 [Google Scholar]
 Benson, D., Stein, R., & Nordlund, Å. 2006, in Solar MHD Theory and Observations: A High Spatial Resolution Perspective, ed. J. Leibacher, R. F. Stein, & H. Uitenbroek, ASP Conf. Ser., 354, 92 [Google Scholar]
 Birch, A. C., & Gizon, L. 2007, Astron. Nachr., 328, 228 [NASA ADS] [CrossRef] [Google Scholar]
 Cameron, R., Gizon, L., & Duvall, Jr., T. L. 2008, Sol. Phys., 251, 291 [NASA ADS] [CrossRef] [Google Scholar]
 ChristensenDalsgaard, J., Dappen, W., Ajukov, S. V., et al. 1996, Science, 272, 1286 [CrossRef] [Google Scholar]
 Couvidat, S., Birch, A. C., & Kosovichev, A. G. 2006, ApJ, 640, 516 [NASA ADS] [CrossRef] [Google Scholar]
 Duvall, Jr., T. L., & Birch, A. C. 2010, ApJ, 725, L47 [NASA ADS] [CrossRef] [Google Scholar]
 Duvall, Jr., T. L., & Gizon, L. 2000, Sol. Phys., 192, 177 [NASA ADS] [CrossRef] [Google Scholar]
 Duvall, Jr., T. L.,Jefferies, S. M., Harvey, J. W., & Pomerantz, M. A. 1993, Nature, 362, 430 [NASA ADS] [CrossRef] [Google Scholar]
 Duvall, Jr., T. L., D’Silva, S., Jefferies, S. M., Harvey, J. W., & Schou, J. 1996, Nature, 379, 235 [NASA ADS] [CrossRef] [Google Scholar]
 Duvall, Jr., T. L.,Kosovichev, A. G., Scherrer, P. H., et al. 1997, Sol. Phys., 170, 63 [NASA ADS] [CrossRef] [Google Scholar]
 Duvall, Jr., T. L.,Birch, A. C., & Gizon, L. 2006, ApJ, 646, 553 [NASA ADS] [CrossRef] [Google Scholar]
 Gizon, L., & Birch, A. C. 2004, ApJ, 614, 472 [NASA ADS] [CrossRef] [Google Scholar]
 Gizon, L., & Birch, A. C. 2005, Liv. Rev. Sol. Phys., 2, 6 [Google Scholar]
 Gizon, L., Duvall, Jr., T. L., & Larsen, R. M. 2000, J. Astrophys. Astron., 21, 339 [Google Scholar]
 Gizon, L., Schunker, H., Baldner, C. S., et al. 2009, Space Sci. Rev., 144, 249 [NASA ADS] [CrossRef] [Google Scholar]
 Gizon, L., Birch, A. C., & Spruit, H. C. 2010, ARA&A, 48, 289 [Google Scholar]
 Hanasoge, S. M., & Duvall, T. L. 2009, ApJ, 693, 1678 [NASA ADS] [CrossRef] [Google Scholar]
 Hansen, P. C. 1998, Rankdeficient and discrete illposed problems: numerical aspects of linear inversion (Philadelphia, PA, USA: Society for Industrial and Applied Mathematics) [Google Scholar]
 Jackiewicz, J., Gizon, L., Birch, A. C., & Thompson, M. J. 2007, Astron. Nachr., 328, 234 [NASA ADS] [CrossRef] [Google Scholar]
 Jackiewicz, J., Gizon, L., & Birch, A. C. 2008, Sol. Phys., 251, 381 [Google Scholar]
 Jackiewicz, J., Birch, A. C., Gizon, L., Hanasoge, S., & Hohage, T. 2011, Sol. Phys., accepted [Google Scholar]
 Jensen, J. M., Jacobsen, B. H., & ChristensenDalsgaard, J. 1998, in Structure and Dynamics of the Interior of the Sun and Sunlike Stars, ed. S. Korzennik, ESA SP, 418, 635 [NASA ADS] [Google Scholar]
 Kosovichev, A. G. 1996, ApJ, 461, L55 [NASA ADS] [CrossRef] [Google Scholar]
 Kosovichev, A. G., & Duvall, Jr., T. L. 1997, in SCORe’96: Solar Convection and Oscillations and their Relationship, ed. F. P. Pijpers, J. ChristensenDalsgaard, & C. S. Rosenthal, Astrophys. Space Sci. Lib., 225, 241 [Google Scholar]
 Moradi, H., Baldner, C., Birch, A. C., et al. 2010, Sol. Phys., 267, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Pijpers, F. P., & Thompson, M. J. 1992, A&A, 262, L33 [NASA ADS] [Google Scholar]
 Rempel, M., Schüssler, M., & Knölker, M. 2009, ApJ, 691, 640 [NASA ADS] [CrossRef] [Google Scholar]
 Scherrer, P. H., Bogart, R. S., Bush, R. I., et al. 1995, Sol. Phys., 162, 129 [Google Scholar]
 Ustyugov, S. D. 2008, in Subsurface and Atmospheric Influences on Solar Activity, ed. R. Howe, R. W. Komm, K. S. Balasubramaniam, & G. J. D. Petrie, ASP Conf. Ser., 383, 43 [Google Scholar]
 Vögler, A., Shelyag, S., Schüssler, M., et al. 2005, A&A, 429, 335 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Zhao, J., & Kosovichev, A. G. 2004, ApJ, 603, 776 [NASA ADS] [CrossRef] [Google Scholar]
 Zhao, J., & Kosovichev, A. G., & Duvall, Jr., T. L. 2001, ApJ, 557, 384 [NASA ADS] [CrossRef] [Google Scholar]
 Zhao, J., Georgobiani, D., Kosovichev, A. G., et al. 2007, ApJ, 659, 848 [NASA ADS] [CrossRef] [Google Scholar]
 Zhao, J., Parchevsky, K. V., Hartlep, T., & Kosovichev, A. G. 2010, in Am. Astron. Soc. Meet. Abstr., 216, 319.05 [Google Scholar]
All Tables
Predicted inversion errors and the expected magnitude of the velocities at the three different depths discussed throughout this paper.
Correlation analysis of various inversion components using synthetic travel times averaged over 4 days.
All Figures
Fig. 1 Examples of synthetic traveltime maps used in this paper (top row) and observed travel times from SOHO/MDI (bottom row). Temporal length of observation is T = 6 h. 

In the text 
Fig. 2 Comparison of azimuthally averaged spatial power spectra of synthetic (solid) and observed (dashed) travel times displayed in Fig. 1 as a function of k R_{⊙}, where k is the horizontal wavenumber and R_{⊙} is the radius of the Sun. 

In the text 
Fig. 3 Example inversion weights for 1 Mm depth and effect of the minimisation of the spread S. A strong regularisation ϵ confines the inversion weights in the spatial domain (bottom row, ϵ = 10^{0.5}) compared to the case ϵ = 0 (top row). The cut at y = 0 (horizontal black line) is displayed on the righthandside panel, where the values of were scaled by . 

In the text 
Fig. 4 The cut through the x = y = 0 point of the averaging kernel (solid) and the respective target function (dashed) for v_{x} inversion with minimised crosstalk, using travel times averaged over 4 days, at three discussed depths (1, 3.5, and 5.5 Mm). 

In the text 
Fig. 5 Comparison of inverted v_{x} with input data. Toprow panels show the input flow field convolved with the target function and the averaging kernel respectively. Bottom rows show inversion results in cases when crosstalk is ignored (left) and minimised (right). Random errors of the inversion are given in Table 2. 

In the text 
Fig. 6 Inversion biases for v_{x} at 1 Mm depth. Left: noiseless with crosstalk minimised versus the ideal . The departure from slope unity (grey dashed line) is due to an imperfect match between and . Right: noisy versus noiseless . The results are plotted in two cases: when the crosstalk is ignored (grey × ) and minimised (black + ). The linear fit to the black crosses coincides with the dashed line of slope unity. The black dotted lines represent the predicted error of 14 m s^{1}, which is consistent with the observed scatter of the black crosses. 

In the text 
Fig. 7 The azimuthallyaveraged power spectra of various components of the v_{x} inversions using travel times averaged over 4 days. For reference, we plot the power spectrum of (thick grey solid line). We plot also the power spectrum of (solid line) and the power spectrum of the noise (i.e., the power spectrum of ; dashed line) for the inversion where the crosstalk is minimised (black) and ignored (grey). 

In the text 
Fig. 8 Horizontal spatial crosscorrelation between v_{x} and v_{z} from the simulation in the top 1 Mm. These flow correlations together with a nonvanishing crosstalk averaging kernel may lead to a systematic bias of the inferred v_{z}. The circle of radius of 15 Mm shows the width of the desired target function. 

In the text 
Fig. 9 The cut through the x = y = 0 point of the averaging kernel (solid) and the respective target function (dashed) for v_{z} inversion using travel times averaged over 4 days. Compare with Fig. 4. 

In the text 
Fig. 10 All components of the averaging kernel for v_{z} inversion at 1 Mm depth with a FWHM of s_{z} = 1.1 Mm and s_{h} = 15 Mm. Bottom row: with crosstalk minimised, top row: crosstalk is ignored. The crosstalk is presented in the form of and averaging kernel components. Random error of the results is 3 m s^{1} when assuming data averaged over 4 days. Overplotted contours, which are also marked on the colour bar for reference, denote the following: halfmaximum of the kernel (white), halfmaximum of the target function (black), and by grey lines ± 5% of the maximum value of the kernel (solid and dotted, respectively). 

In the text 
Fig. 11 Comparison of inverted v_{z} with input data, compare with Fig. 5. ⟨ v_{z} ⟩ = 11 m s^{1}. 

In the text 
Fig. 12 All components of the v_{z} inversion at 1 Mm depth. Top row with the crosstalk ignored, bottom row with the crosstalk minimised. We demonstrate that if crosstalk is not addressed, horizontal components will leak into the inverted v_{z} and cause a bias. 

In the text 
Fig. 13 Inversion biases for v_{z} at 1 Mm depth, similar to Fig. 6. While the black line (with minimised crosstalk) almost coincides with the dashed line of slope unity, the grey one (crosstalk ignored) indicates a horrible bias. The predicted error is 3 m s^{1}. The inversion is not sensitive to the horizontally averaged horizontal flow ⟨v_{z}⟩ = 11 m s^{1}, therefore it was subtracted from in this comparison. 

In the text 
Fig. 14 The azimuthallyaveraged power spectra of various components of the v_{z} inversions using travel times averaged over 4 days. Compare with Fig. 7 for v_{x} inversions. Notice that in the case of containing the crosstalk at 1 Mm depth, the excess in power around kR_{⊙} ~ 150 is not caused by the random noise, but by the bias coming from the crosstalk (not plotted separately). In the other two cases the random noise is a main cause why these inversions are not possible. 

In the text 
Fig. 15 Comparison at three depths of statistical averages of velocityvector over many flow realisations. The horizontal flow is displayed by arrows while the vertical flow is colourcoded. The reference arrow indicates 100 m s^{1}, the random errors of the inversions are given in Table 2. Compare to Figs. 5 and 11. 

In the text 
Fig. 16 Slices through the x and z components of the simulated flow v at depths 1 Mm, 3.5 Mm, and 5.5 Mm. 

In the text 
Fig. 17 The pointtoannulus sensitivity kernels for three flow components computed for the fmode, distance 28 Mm and outwardinward geometry. The white circle represents the location of the averaging annulus. 

In the text 
Fig. 18 Horizontal averages of sensitivity kernels may serve useful when estimating which depths are easier to target. The trend for each mode/ridge was obtained by taking and averaging over all as within the given mode. 

In the text 
Fig. 19 An example noisecovariance matrix for fmode travel times averaged over 6 h. In this plot, a stands for the combination of the fmode, oi geometry, and annulus radius of 7.3 Mm, b stands for the combination of fmode, we geometry, and annulus radius of 8.8 Mm. 

In the text 
Fig. 20 All components of the averaging kernel for v_{x} inversion at 1 Mm depth with a FWHM of s_{z} = 1.1 Mm and s_{h} = 15 Mm. Bottom row: with crosstalk minimised, top row: crosstalk is ignored. Random error of the results is 14 m s^{1} when assuming data averaged over 4 days. Overplotted contours, which are also marked on the colour bar for reference, denote the following: halfmaximum of the kernel (white), halfmaximum of the target function (red), and ± 5% of the maximum value of the kernel (blue and green, respectively). 

In the text 
Fig. 21 All components of the averaging kernel for v_{x} inversion at 3.5 Mm depth with a FWHM of s_{z} = 2.2 Mm and s_{h} = 15 Mm. Random error of the results is 20 m s^{1} when assuming data averaged over 4 days. For details see Fig. 20. 

In the text 
Fig. 22 All components of the averaging kernel for v_{x} inversion at 5.5 Mm depth with a FWHM of s_{z} = 3.5 Mm and s_{h} = 15 Mm. Random error of the results is 28 m s^{1} when assuming data averaged over 4 days. For details see Fig. 20. 

In the text 
Fig. 23 The contributions of particular modes to the horizontally averaged averaging kernel for v_{x} inversions using travel times averaged over 4 days for depths 1 and 3.5 Mm. We do not display the inversion for the depth of 5.5 Mm, because it is heavily dominated by noise. 

In the text 
Fig. 24 To solve the peculiarity of the v_{z} inversion, we introduced two formalisms in Sect. 5.1. Here we plot performance of those. In black, the magnified section a y = 0 and z = z_{0} of different target functions are displayed, one with removed mean (1) and one constructed with negative sidelobes (2). The resulting averaging kernels are also plotted. It is evident that the resulting averaging kernels are qualitatively very similar even when different formalisms were used to compute them. 

In the text 
Fig. 25 All components of the averaging kernel for v_{z} inversion at 3.5 Mm depth with a FWHM of s_{z} = 2.2 Mm and s_{h} = 15 Mm. Random error of the results is 13 m s^{1} when assuming data averaged over 4 days. For details see Fig. 20. 

In the text 
Fig. 26 All components of the averaging kernel for v_{z} inversion at 5.5 Mm depth with a FWHM of s_{z} = 3.5 Mm and s_{h} = 15 Mm. Random error of the results is 133 m s^{1} when assuming data averaged over 4 days. For details see Fig. 20. 

In the text 
Fig. 27 The cut through the x = y = 0 point of the averaging kernel (solid) and the respective target function (dashed) for the v_{x} (left) and v_{z} (right) inversions using averaging over many flow realisations plotted along with the corresponding target functions at three discussed depths (1 Mm in blue, 3.5 Mm in green, and 5.5 Mm in red). Compare to Figs. 4 and 9 where the resemblance of the target functions is worse. The random error of the results is given in Table 2. 

In the text 
Fig. 28 The azimuthallyaveraged power spectra of the v_{z} inversion components at depths of 1, 3.5 and 5.5 Mm for averaging over many flow representations. For reference, we plot the power spectrum of using the black solid line. Then we plot the power spectrum of (solid line) and power spectrum of the noise (i.e., the power spectrum of ; dashed line) for the inversion where the crosstalk is minimised (blue) and ignored (red). Compare to Fig. 14. 

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.