Issue 
A&A
Volume 622, February 2019



Article Number  A163  
Number of page(s)  12  
Section  The Sun  
DOI  https://doi.org/10.1051/00046361/201833000  
Published online  14 February 2019 
Combined helioseismic inversions for 3D vector flows and soundspeed perturbations
^{1}
Astronomical Institute of Charles University, Faculty of Mathematics and Physics, V Holešovičkách 2, Praha 8 180 00, Czech Republic
email: d.korda@seznam.cz
^{2}
Astronomical Institute of Czech Academy of Sciences, Fričova 298, Ondřejov 25165, Czech Republic
Received:
10
March
2018
Accepted:
20
December
2018
Context. Time–distance helioseismology is the method of the study of the propagation of waves through the solar interior via the travel times of those waves. The travel times of wave packets contain information about the conditions in the interior integrated along the propagation path of the wave. The travel times are sensitive to perturbations of a variety of quantities. The usual task is to invert for the vector of plasma flows or the sound–speed perturbations separately. The separate inversions may be polluted by systematic bias, for instance, originating in the leakage of vector flows into the sound–speed perturbations and vice versa (called a crosstalk). Information about the crosstalk is necessary for a proper interpretation of results.
Aims. We introduce an improved methodology of the timedistance helioseismology which allows us to invert for a full 3D vector of plasma flows and the sound–speed perturbations at once. Using this methodology one can also derive the mean value of the vertical component of plasma flows and the crosstalk between the plasma flows and the sound–speed perturbations.
Methods. We used the Subtractive Optimally Localised Averaging method with a minimisation of the crosstalk as a tool for inverse modelling. In the forward model, we use Born approximation traveltime sensitivity kernels with the Model S as a background. The methodology was validated using forwardmodelled travel times with both mean and difference pointtoannulus averaging geometries applied to a snapshot of fully selfconsistent simulation of the convection.
Results. We tested the methodology on synthetic data. We demonstrate that we are able to recover flows and sound–speed perturbations in the nearsurface layers. We have taken the advantage of the sensitivity of our methodology to entire vertical velocity, and not only to its variations as in other available methodologies. The crosstalk from both the vertical flow component and the sound–speed perturbation has only a negligible effect for inversions for the horizontal flow components. Furthermore, this crosstalk can be minimised if needed. The inversions for the vertical component of the vector flows or for the sound–speed perturbations are affected by the crosstalk from the horizontal components, which needs to be minimised in order to provide valid results. It seems that there is a nearly constant crosstalk between the vertical component of the vector flows and the sound–speed perturbations.
Key words: Sun: helioseismology / Sun: oscillations / Sun: interior
© ESO 2019
1. Introduction
The internal structure of the Sun and the dynamics of plasmas therein are key factors in understanding the solar activity, which in turn influences the interplanetary space (including effects on human infrastructures). Even though processes in the interior of the Sun may be modelled using numerical codes (e.g. Brun et al. 2004; Rempel et al. 2009), these models must be constrained by observations to ensure that they describe a physical reality to some extent.
Helioseismology is the only method that allows us to study subsurface layers of the Sun; these layers, at least to some extent, reflect mostly the properties of the convection zone. Knowledge of the structure and dynamics of plasmas in these regions allows us to put important constraints on theories of the solar dynamo and of the formation and evolution of magnetic fields, including the models of sunspots, and others. The models and inferences from helioseismology do not always agree. For instance, the estimate of convective velocity in the convection zone from helioseismology is two orders of magnitude less than that derived from theoretical studies (Hanasoge et al. 2012). These results, however, were not confirmed by an independent investigation by Greer et al. (2015) when using a different helioseismic analysis. Also, from theoretical considerations (Miesch et al. 2012) the results of Hanasoge et al. (2012) do not fit into the current view of convection in the Sun. On the other hand, these discrepancies have not been explained until now. Another surprising result is that there are helioseismic indications that the supergranules may have a different depth structure than predicted by theory (Duvall & Hanasoge 2013; Švanda 2012; Duvall et al. 2014). Even though both mentioned results are considered to be controversial, there is an obvious need for both theory and helioseismology to reconcile these existing discrepancies.
Helioseismology is used to study the solar interior via interpretation of the dispersion relations of acoustic and surface gravity waves. These waves are generated randomly via the vigorous convection in the convection zone. The waves manifest themselves in many different observables, such as intensity maps of the photosphere or Doppler shifts of photospheric spectral lines. The focus of global helioseismology is to infer the background structure of the solar interior. Local helioseismology focuses on the detection and quantification of small localised perturbations of plasma parameters with respect to the reference solar model.
In this study, we focused on the second approach, using the time–distance method (Duvall et al. 1993) to study the perturbations of travel times of waves caused by inhomogeneities of plasma parameters about the reference model. The time–distance method has been a standard method of solar research for several decades. A comprehensive review of methods of local helioseismology, including the time–distance method and results obtained by this method was presented for example in Gizon et al. (2010) or in Kosovichev & Zhao (2016). The time–distance helioseismology has been used in many applications: for inversions for the solar rotation (Schou 1997), meridional circulation (Zhao & Kosovichev 2004), perturbations of density, pressure, temperature, sound speed (Tong et al. 2003; Brüggen & Spruit 2000), and plasma flows (e.g. Švanda 2013).
The usual task of time–distance helioseismology is to invert for the horizontal components and/or spatial variations of the vertical component of the vector of plasma flows, or the sound–speed perturbations separately, using different sets of traveltime measurements. The separate inversions are justified under the assumption of independent action of various perturbers to the wave travel times. In the inverse modelling, however, it is not certain that this assumption holds. Also, in the realistic models of solar plasmas, the perturbers are correlated. For instance, in the massconserving convective flows, the upflows (in the vertical direction) in the cell centres are accompanied by outflows (in the horizontal direction). Thus, these two quantities are naturally correlated. This correlation between horizontal and vertical velocities was convincingly shown by Švanda et al. (2011); see their Fig. 8. A similar correlation plot between the vertical velocity and the sound–speed perturbation is given in Fig. 1, where a large correlation (correlation coefficient of 0.48 at 0.5 Mm depth) between these two quantities is seen. The correlation coefficient was calculated from a realistic solar convection simulation provided by Rempel (2014). This correlation has its origin in the fact that the large upflow velocities are usually located in the interior of the convective cells, where the temperature is greater compared to the surroundings. The greater temperature implies an increase in the sound speed, which naturally follows from the idealgas equation of state.
Fig. 1.
Map of correlation between the vertical velocity and the sound–speed perturbation in the nearsurface layers obtained from the realistic simulation of the solar convection by Rempel (2014). This plot is analogous to Fig. 8 of Švanda et al. (2011). The black circle corresponds to the width of the target function (15 Mm). 

Open with DEXTER 
In time–distance helioseismology, the travel times are usually averaged over the annulus around the cell centre, which, together with the correlation described above, leads to a leakage – the crosstalk – of one perturber into the inversion for another one. Separate inversions do not allow one to quantify the possible crosstalk contributions. Information about the crosstalk is extremely important for proper interpretation of results, especially for the vertical velocity (Švanda et al. 2011). This is because the amplitudes of variations of the vertical velocity are less than 10 m s^{−1}, but the amplitudes of the horizontal velocities are on the order of 100 m s^{−1}. Even a small crosstalk can absolutely devalue the results.
Our aim is to combine inversion for the vector flows υ = (υ_{x}, υ_{y}, υ_{z}), where (x, y, z) are local Cartesian coordinates, and for the sound–speed perturbations δc_{s}. This approach will allow us to study, for the first time, the crosstalk between the flow perturbations and the sound–speed perturbations. The expected amplitude of the sound–speed perturbations known from the numeric simulations is approximately 50 m s^{−1}, which is comparable with the amplitude of the vertical velocity, where the crosstalk has a very important effect.
Another remaining issue in the description of solar flows is that with the current methodology it is possible to invert only for the spatial variations of the vertical velocity. Due to the symmetries of the sensitivity kernels used in the flow inversions, the mean value of the vertical component vanishes. The improved methodology we propose in this paper, that is, merging inversions for the 3D flows and the sound–speed perturbations into one inversion, allows us to fully describe the vertical flow component.
2. Methodology
2.1. Forward problem
We measure the perturbed travel time δτ by crosscorrelating the signal (the Doppler shift of a photospheric spectral line) in a point at a horizontal location r with the signal averaged over the surrounding annulus with the radius Δ; we use the approach of Gizon & Birch (2004). For each Δ, four different geometries are introduced by different weightings in the azimuthal angle (three difference types and one mean type). We use the outflow–inflow geometry denoted by o–i, which is the difference between the travel times of the waves travelling from the annulus centre to its rim and the waves travelling in the opposite direction. The east–west and north–south, denoted by e–w and n–s, are computed similarly, only the signal over the annulus is weighted by either cosine (e–w) or sine (n–s) of the polar angle about the central point. Thus, the measurement is sensitive to the waves travelling in one of these directions. The last one is the mean geometry, which averages travel times of the waves travelling from the centre to the rim of the annulus and of the waves travelling in the opposite direction. In the inversions seen in the literature, the difference geometries were used only to invert for flows, whereas the mean geometry was used to invert for the sound speed.
The perturbed travel time δτ can theoretically be computed from the model as
where r and r′ are the horizontal positions, z is the vertical position, and are the sensitivity kernels from the forward modelling (Gizon & Birch 2004; Birch & Gizon 2007). Here we use Model S (ChristensenDalsgaard et al. 1996) as a reference background model. The δq_{β} are the perturbations of the physical quantities (such as flows, sound–speed, density, etc.; indexed by β), P is the number of these quantities, and n^{a} is a realisation of the random noise (Jackiewicz et al. 2012). The superscript a corresponds to the measurements of individual travel times and associates the selection of wave filters and averaging geometries as described above, including a selection of the radius of the annulus.
2.2. Inverse problem
The goal of an inverse problem is to invert Eq. (1) and compute
at a position (r_{0};z_{0}), where α indicates the perturber (for instance the vertical flow velocity). For simplicity we use a subscript α to indicate the quantity in the direction of the inversion and a subscript β as an index indicating a general physical quantity component. The above illustrated inversion is generally not possible because of high levels of noise. The time–distance helioseismology methods are used to derive an estimate of the given quantity, denoted as .
In our methodology we focus on a method called Optimally Localised Averaging (OLA), which corresponds to the method of an approximate inverse in mathematical terminology (Louis & Maass 1990). The OLA method (Backus & Gilbert 1968) was originally developed for Earth seismology but can also be used for helioseismology (Pijpers & Thompson 1992).
The SOLA method. For our study we use an improved version of the Subtractive OLA method (SOLA) described by Jackiewicz et al. (2012). The SOLA searches for at a given position (r_{0};z_{0}) in the form of a linear combination of the travel times and the unknown weight functions w^{α}:
where N is the total number of the horizontal positions and M is the number of the traveltime geometries. The functions w^{α} minimise the χ^{2} in the form
We use , where 𝒯 is a userselected target function in the direction of the inversion and is a Kronecker delta. The target function is localised around the area of investigation. Λ^{ab}(r_{i}−r_{j}) = cov[n^{a}(r_{i}),n^{b}(r_{j})] describes a noise covariance matrix between observations denoted by a and b in positions denoted by r_{i} and r_{j}. Parameters μ, ν, and ϵ are userselected tradeoff parameters. The quantity is defined as
and is referred to as the averaging kernel. It quantifies the level of smearing of the real quantity δq_{β}. By using the averaging kernels, the real quantity and the estimate from the inversion are linked by
The first term of Eq. (4) is the misfit term between the userselected localisation of the inverted quantity and the real localisation given by . The second term is the regularisation of the propagation of the random noise. The third term is the regularisation of the crosstalk between the inverted quantity and the other quantities and its meaning is obvious from Eq. (6), where the quantity in the direction of the inversion and also other quantities contribute to the inverted estimate. The fourth term minimises the spatial spread of the weight functions and the last term is the normalisation of the averaging kernels added by the Lagrange multipliers λ^{β} (for details see Švanda et al. 2011).
The minimisation of Eq. (4) with respect to w^{α} and λ^{β} gives the solution for the weights. We assume the spatial homogeneity of the background model and therefore the problem is solved in the Fourier domain. There, it decouples and allows us to split a large realistic realspace problem into a series of small Fourierspace problems. This approach permits a realistic inverse problem to be solved even with a desktop computer. One can find the Fourierspace formulation in Jackiewicz et al. (2012) or Švanda et al. (2011).
Noise matrices. We use datadriven estimates for the noise matrices Λ^{ab}. They were computed by using a 1/T fitting introduced by Gizon & Birch (2004), where T is the traveltime averaging time. We used a large set of traveltime maps measured for 31 days by the HMI/SDO instrument (Scherrer et al. 2012; Schou et al. 2012) in the quiet Sun period from August 16, 2010, to September 15, 2010. Following Fournier et al. (2014) we fitted covariances of the travel times δτ^{a} and δτ^{b} by using the function cov[δτ^{a},δτ^{b}] = c + d/T for T = {2, 3, 4, 6, 8, 12, 24} hours, where c and d are the fitted parameters. The parameter d is then the element of Λ^{ab}. We also tried to fit the functions d/T + e/T^{2} and c + d/T + e/T^{2} with the fitting parameters c, d, and e as suggested by Fournier et al. (2014), but the results were not as expected. For instance, the diagonal terms Λ^{aa} (the RMS of the travel times) were overly smeared and did not correspond to the measured RMS values.
The full covariance matrix is important for getting a reasonable estimate of the variation of the inverted quantity , which equals
3. Improvements
The methodology described above is universal and has been used in various inversions. For instance, it was used to study the vector flows (Švanda et al. 2011, and following works) and was found useful for sound–speed inversions (Jackiewicz et al. 2012). These two types of inversions differ in the selection of the appropriate sensitivity kernels in the inverse problem and hence the corresponding traveltime measurements. For the flows the difference geometries were used (i.e. the o–i, e–w, and n–s geometries), whereas for the sound–speed inversions the mean geometry was used. This means that the subscript β in Eq. (1) incorporated either the vector of the flow only with δq_{β} = (υ_{x},υ_{y},υ_{z}) or the fractional sound–speed perturbations only with δq_{β} = ().
This is in agreement with recent findings by Burston et al. (2015) where the authors show in their Table 1 that the horizontal flow components υ_{x} and υ_{y} are largely sensitive to the difference types of the travel times (o–i, e–w, and n–s), whereas δc_{s} is sensitive to the mean type. The vertical flow component υ_{z} is sensitive to both types however; to invert for the full υ_{z} and not only its spatial variations, the meantype geometry needs to be included.
Inspired by such findings, we modified the working methodology so that the pipeline can be used to perform inversions for the vector flows and the sound–speed perturbations at once. In effect, the subscript β in Eq. (1) now enumerates all components of the flow and the sound–speed perturbations at once, β = (x,y,z,s). Thus, the vector δq_{β} = (υ_{x},υ_{y},υ_{z},δc_{s}).
Consequently, the traveltime geometries contain all four averaging geometries, the radii of the averaging annuli, and mode filtering – all of these choices are covered in the aggregate superscript a in the equations.
This approach allows us to invert for the complete vector of the plasma flows and the sound–speed perturbations in one inversion, which has not been done before.
We note that variables were selected after consideration of extensive testing of four possible constructions of the vector δq_{β}:

1.
The first one was (υ/c_{s},δc_{s}/c_{s}). There is a problem with the interpretation of results, because in order to convert the results into physical units, one has to multiply them by the background sound–speed profile, smeared by the averaging kernel. The averaging kernels for the inverted horizontal velocities scaled by the sound speed had a deeper secondary lobe that could not be removed by varying the tradeoff parameters, which made the interpretation of results difficult. Because of the deeper lobe the horizontal velocities were not localised around the target depth. The noise estimates targeted by the inversion were biased.

2.
Another option was (υ,δc_{s}), which we took as the preferred one. This option has many advantages: The quantities are in physical units. The maps of the quantities, the crosstalk, and the localisation of are easy to interpret. The averaging kernels do not have the secondary lobe like in the previous case and are welllocalised around the target depth.

3,4.
The other options were (υ,/c_{s}) and (υ/c_{s},). If we selected one of these options we would have to convert /c_{s} or υ/c_{s} and to the physical units similarly to option #1.
The traveltime kernel codes^{1} usually compute the kernels for the fractional squared perturbation of the sound speed, that is, for . Instead of modifying and extensively testing the kernel code, we rather compute the kernels for linear δc_{s} from the fractional squared perturbations using the following approach.
We write the squared sound–speed perturbation as
where c_{model} is the background sound–speed value. The second term has a negligible magnitude compared to the first term. We then compute δc_{s} sensitivity kernels in the form
where we neglected against (the uncertainty is less than 1%). Finally, we obtain the relation for the sensitivity kernel for δc_{s}:
4. Validation using synthetic data
4.1. Synthetic travel times
We validate our new methodology using synthetic data. A snapshot from a box^{2} of realistic simulation of the convection zone (Rempel 2014; DeGrave et al. 2014) was used. It naturally contains both the complete vector of the plasma flows and the sound–speed perturbations.
The synthetic forwardmodelled travel times were computed using Eq. (1). We further add a realistic representation of the random noise to the travel times. Having these two separate components, we were able to study various signaltonoiseratio (S/N) situations.
The realisation of the random noise was computed in the Fourier space using a generator of multivariate normal random numbers using the covariance for each wave vector k separately, where is the Fourier transform of Λ^{ab}. We forced the expected Fourierspace symmetries so that in the real space the noise realisation had only a real component as expected. By running a large number (10 000) of realisations and computing the noise covariance matrix again we verified that the random noise realisations have the same statistical distribution as the corresponding Λ^{ab}.
4.2. Comparison with the model
We applied our improved pipeline to the synthetic travel times and compared the results to the known inputs. By using the synthetic data, we have a complete understanding of the contributions of the individual components of the inversions, that is, all components of the averaging kernels, the noise, and also the expected ideal answer. Thanks to the use of the synthetic data, all these components may be validated separately for each quantity δq_{α} = {υ_{x}, υ_{y}, υ_{z}, δc_{s}}. The inversion components to be validated are
The ideal answer is what we expect to obtain by applying the inversion. The deviations between the ideal answer and our result are caused by the random noise and the otherthanexpected localisation, described by the 4D averaging kernel.
For this particular study, we used only the fmode ridgefiltered travel times. The reason for this is the strong localisation towards the surface of the Sun, which means that the results, if applied to the real data, can directly be confronted with the results obtained by other methods (similarly to e.g. Švanda et al. 2013). The disadvantage of using the f mode is that it is not particularly sensitive to the sound–speed perturbations (Burston et al. 2015). According to the same paper, the acoustic p modes should have an orderofmagnitude better sensitivity to the sound–speed perturbations. To confirm this we performed equivalent inversions using the p_{1} mode. These results are discussed in Sect. 7 (Validation with the p_{1} mode).
4.3. The target function
In the SOLA inversions, the user is free to select the target function 𝒯, which represents the firstorder estimate of the averaging kernel. We describe the results of the comparison with a GaussiandampedGaussiantype (GdG) target function, which consists of the 3D Gaussian with the horizontal full width at half maximum (FWHM) of 15 Mm and the vertical FWHM of 0.5 Mm. The centre of the target function is placed at a depth of 0.5 Mm. In the vertical direction, the target function is further multiplied by a smooth function so that 𝒯 reaches zero at the solar surface, z = 0, and does not extend to the domain of the solar atmosphere.
4.4. Tradeoff parameters
The inverted depends on a selection of the tradeoff parameters μ, ν, and ϵ. We selected the tradeoff parameters in the following manner: At first we set an upper limit of the level of the random noise, according to the inverted quantity (20 m s^{−1} for υ_{x} and υ_{y}, 15 m s^{−1} for δc_{s}, 5 m s^{−1} for υ_{z}). This gave us a lower limit of the parameter μ.
In the second step we plotted the averaging kernels and cuts through the weights for all combinations of the tradeoff parameters which fulfilled the noise constraint. Our further aims were to minimise the crosstalk which is regularised by the parameter ν. The parameter ϵ controls the spatial spread of weights to ensure that the inversion weights are spatially confined and to avoid ringing solutions. Larger values of the tradeoff parameters are usually balanced by a worse fit of the target function by the averaging kernel. The optimal combination of tradeoff parameters is chosen on the trialanderror basis. We note that there is an Lcurve (e.g. Hansen 1999) method which attempts to determine the values of the tradeoff parameters in an exact way. Unfortunately, in our case this method usually leads to a large level of random noise and therefore we are not using it.
4.5. The method of evaluation
In the following sections we discuss the inversions for the horizontal and vertical flow and the sound–speed perturbations separately. When applicable, we compare the new results to the outputs of the Švanda et al. (2011) pipeline.
For all perturbations we always discuss the shape of the resulting averaging kernel to estimate the localisation in the Sun, the smearing, and the crosstalk contributions. In the plots we show the horizontal slices through the target depth and also the vertical slices perpendicular to the symmetry axis. We then describe the individual contributions to the inverted quantity and discuss the fulfilment of the inversion requirements, such as the spatial localisation of the weights.
5. Inversions for horizontal components of vector flows
Inversions for horizontal (that is υ_{x} and υ_{y}) components of the vector flows have been typical tasks for the timedistance helioseismology. There is only a small amount of crosstalk with the vertical υ_{z} caused by natural correlations between the vertical upflow and the outflows in the massconserving flow. The averaging kernels are plotted in Fig. 2 (crosstalk minimised) and Fig. 3 (crosstalk not minimised). We note that only plots for υ_{x} are given, plots for υ_{y} are similar except for the rotation by 90 degrees around the vertical axis. By comparing the second, third, and fourth columns of the two figures, one can see that the crosstalk minimisation works as expected, when these contributions vanish in Fig. 2. The crosstalk minimisation is balanced by the worse fit of the component of the averaging kernel in the direction of the inversion to the target function. The averaging kernel is more extended in the x direction, reaching essentially an elliptical shape in the horizontal cut, whereas the target function is roundish. The use of the fmode sensitivity kernels only also contributed to this issue because they do not represent the GdG function well.
Fig. 2.
Averaging kernels for υ_{x} inversion at 0.5 Mm depth. The crosstalk was minimised. The red curve corresponds to the halfmaximum of the GdG target function, the yellow curve corresponds to the halfmaximum of the averaging kernel and the blue solid and blue dotted (not present in this figure) curves correspond to +5% and −5% of the maximum of the averaging kernel, respectively. In the top row there are horizontal slices of the averaging kernel at the target depth and in the bottom row there are vertical slices perpendicular to the symmetries. In the right part one can see the selection of the tradeoff parameters and the standard deviation when assuming the travel times averaged over 24 hours. 

Open with DEXTER 
Fig. 3.
Similar to Fig. 2: averaging kernels for υ_{x} inversion at 0.5 Mm depth without crosstalk minimisation. 

Open with DEXTER 
In Fig. 4 one can see the individual contributions from the different quantities and the noise to both inversions (with and without the minimisation of the crosstalk). There is also a comparison with the ideal answer. Despite the appearance of the averaging kernel with crosstalk not being minimised, the total leakage from υ_{z} and δc_{s} is negligible because the crosstalk components of the averaging kernel have a rather small extent. Moreover, the amplitudes of the leaking components are at least an order of magnitude smaller in the nearsurface layers.
Fig. 4.
Top row: inversion with minimisation of the crosstalk. Bottom row: inversion without minimisation of the crosstalk. First column: contribution to . Second column: contribution to . Third column: contribution to . Fourth column: contribution to . Fifth column: noise contribution to . Sixth column: . Seventh column: ideal answer. 

Open with DEXTER 
Due to the fact that the amplitudes of the horizontal flows are at least an order of magnitude larger than the amplitudes of the two remaining perturbers, there is no need to use the improved methodology for such inversions. It might be unnecessary to minimise the crosstalk at all. Such a claim will, however, lose its validity in the layers or regions, where the horizontal and vertical components reach a comparable magnitude.
6. Inversion for the vertical flow component
Until now, only variations of the vertical component υ_{z} were known. This is because only the difference traveltime geometries were used to invert for υ_{z} and the total horizontal integral of the sensitivity kernels vanished due to their symmetries (Burston et al. 2015). Therefore, such measurements are not sensitive to the horizontal average (the mean) ⟨υ_{z}⟩ which is therefore identically zero in the inverted estimates. We only note that this property also requires a modification of the target function for such an inversion (Švanda et al. 2011) in order to fulfil the normalisation constraint in Eq. (4).
According to Burston et al. (2015) the mean traveltime geometry is sensitive to ⟨υ_{z}⟩. Therefore, by adding the mean geometry we can invert for full υ_{z}. From the previous works (e.g. Švanda et al. 2011) it is known that there is a large crosstalk with the horizontal flow components, which must be minimised in order to have credible results.
The inversion averaging kernels are plotted in Figs. 5 (crosstalk minimised) and 6 (crosstalk unconstrained). The crosstalk minimisation approach allows us to minimise the leakage from the horizontal flow components, however, it has only a limited effect on the leakage of the sound–speed perturbations. We found that this crosstalk with δc_{s} has the origin in the contribution of the mean traveltime kernels. To estimate the severity of this crosstalk, we plotted the meangeometry contribution of the four physical quantities to the averaging kernel. As one can see from Fig. 7, the mean contribution comes from the kernel, with a small contribution from the kernel. There is no contribution from the horizontal flow, because the mean traveltime sensitivity kernels for the horizontal components have vanishing horizontal integrals. From this plot we estimate that the leakage from the sound–speed perturbation to the inversion for the vertical flow is very small.
Fig. 5.
Averaging kernels for υ_{z} inversion at 0.5 Mm depth. The crosstalk was minimised. See Fig. 2 for details. 

Open with DEXTER 
Fig. 6.
Averaging kernels for υ_{z} inversion at 0.5 Mm depth. No crosstalk minimisation. See Fig. 2 for details. 

Open with DEXTER 
Fig. 7.
Horizontal integral of the mean geometry contribution of various physics quantities to the averaging kernel. 

Open with DEXTER 
The effect of the minimisation of the crosstalk with the horizontal component is clearly demonstrated in Fig. 8, where individual inversion components are plotted. The crosstalk contributions seen in the bottom row have a significant impact on the sum and their value in the case of the horizontalcomponents leakage is higher than the typical variations of . For a better comparison between the components and the ideal answer, in Fig. 9 we plot with and without the minimisation of the crosstalk, , and the ideal answer. These plots show that the crosstalk minimisation indeed improves the inversion results. The correlation coefficients between the inverted value and the ideal answer is 0.36 when the crosstalk is ignored and 0.57 when it is minimised. By comparing our new results to the results of Švanda et al. (2011) we can say that even the supplement of the mean travel times in the vertical flow inversion helped to lower the crosstalk, because without them, the crosstalk unconstrained inversions ended with a correlation coefficient of close to −1 with the expected answer (also Zhao et al. 2007).
Fig. 8.
Individual contributions into and its comparison with the ideal answer. See Fig. 4 for details. 

Open with DEXTER 
Fig. 9.
Top row: with minimisation of the crosstalk and without minimisation of the crosstalk. Bottom row: and the ideal answer. The correlation coefficient between and the ideal answer is 0.87. The colour bars are not the same. 

Open with DEXTER 
Obviously, the inversion scales the values by a factor of 0.44, which is due to a large misfit between the target function and the averaging kernel. Such an issue may be understood as being due to the fact that the inferred vertical flows are smoother than those expected from the setup. Such a result highlights the importance of having the information about the averaging kernel when interpreting the inversions.
The incorporation of the mean travel times allowed us to invert for the full vertical velocity and not only its horizontal variations. To demonstrate this, we additionally inverted for the vertical velocity using the old pipeline of Švanda et al. (2011) using only the difference traveltime kernels. It is necessary to modify the target function so that it has a vanishing horizontal mean, which is achieved by introducing a wide negative horizontal side lobe. The results shown in Fig. 10 demonstrate that the expected ideal answers obtained by using two different target functions are essentially indistinguishable, except for the systematic offset. By subtracting the two maps (seen in the third row of Fig. 10) we see that this offset actually has a structure with a very small variation compared to the actual variations of the ideal answer. The offset represents a largescale structure of the vertical flow, which is not detectable using the Švanda et al. (2011) approach.
Fig. 10.
Top row: results of the inversion for υ_{z} with the mean and the difference geometries. Middle row: results of the inversion with the difference geometries only. Bottom row: difference of two aforementioned inversions. Left column: . Right column: ideal answers of the inversions. We point out that the colour bars are not scaled. 

Open with DEXTER 
When the mean travel times are incorporated into the inversion, the overall structure of the inferred vertical flow does not change. The correlation coefficient of the improved vertical flow inversion and the inversion by the original pipeline is 0.91. However the largescale offset is successfully measured by the improved methodology as one can see from the differences in the bottom row of Fig. 10. Except for the lowered amplitude (due to the imperfect averaging kernel), its structure is measured reasonably well.
7. Inversion for soundspeed perturbations
An example of the inversion for the sound–speed perturbations using principally the same methodology we use in this study was published by Jackiewicz et al. (2012). Their inversion used only the mean traveltime sensitivity kernels and thus did not take into account possible leakage from other quantities. We have the opportunity to do this now.
The sound–speed inversion averaging kernels are plotted in Figs. 11 (with the crosstalk minimisation) and 12 (without the crosstalk minimisation). The fits in the direction of the inversion (the rightmost panels) are not very representative of the imposed GdG target function, which is the consequence of using the f mode only. On the other hand, the crosstalk minimisation performs well by essentially setting the crosstalk components of the averaging kernels to zero. For that one has to expect a worse fit of the target function in the direction of the inversion, where even an extended smallamplitude surrounding annulus appears. We assume that the fit improves when more independent traveltime measurements are used in the inversion because of a larger number of independent observations. A wider averaging kernel causes a largerthanexpected smoothing of the results, which is demonstrated by the smaller amplitude of the inverted sound–speed perturbations compared to the ideal answer.
Fig. 11.
Averaging kernels for δc_{s} inversion at 0.5 Mm depth. The crosstalk was minimised. See Fig. 2 for details. 

Open with DEXTER 
Fig. 12.
Averaging kernels for δc_{s} inversion at 0.5 Mm depth. No crosstalk minimisation. See Fig. 2 for details. 

Open with DEXTER 
In Fig. 13 we plot individual contributions to . In both cases the inverted sound–speedperturbations maps are dominated by the component in the direction of the inversion, which has the RMS value of 31 m s^{−1}. The RMS of the total crosstalk is 16 m s^{−1} when not minimised. The crosstalk reaches half of the magnitude of the desired signal, which we consider too large. After the minimisation, the crosstalk contribution lowers to only less than 4 m s^{−1}. The crosstalk is positively correlated with the inverted value, that is, without its minimisation, we overestimate the real sound–speed perturbations by a factor of about 1.9.
Fig. 13.
Individual contributions from and noise to and its comparison with the ideal answer. See Fig. 4 for details. 

Open with DEXTER 
In Fig. 14 one can compare with and without the minimisation of the crosstalk, , and the ideal answer. The correlation coefficient between and the ideal answer is 0.90. The difference between (top left) and the ideal answer is mainly caused by the misfit because the sensitivity of the f mode to δc_{s} is weak.
Fig. 14.
Top row: with minimisation of the crosstalk and without minimisation of the crosstalk. Bottom row: and ideal answer. The correlation coefficient between and ideal answer is 0.90. The colour bars are not the same. 

Open with DEXTER 
7.1. Validation with the p_{1} mode.
According to Burston et al. (2015), the sensitivity of the acoustic modes is at least one order of magnitude higher than that of the f mode. Therefore we ran the sound–speed inversion again, using only the kernels for the p_{1} ridge. The first acoustic ridge samples similar depths to the f mode, and therefore the inversions are directly comparable.
Furthermore, in the inversions presented above, we imposed a certain shape (the “GdG”) of the target function. We showed that in some cases the inverted results were far from optimal, which was due to the fact that the fmode traveltime sensitivity kernels were not dependent enough on depth to fit the GdG target function. For this reason we recomputed the inversion with a target function that has a more natural profile in the vertical domain.
In our methodology, such a target function (called Arbitraryz, Az) resembles the natural depth profile of the sensitivity kernels used in the inversion by representing their kineticenergy profile. It is usually obtained by simply averaging the horizontal integrals as a function of depth of all sensitivity kernels used in the given inversion. Using the Azclass target function our inversion is almost equivalent to the 2D inversion and thus it is expected that the resulting averaging kernels will better fit this target function, while having a lower level of random noise at the same time. In the horizontal direction the target function is again a Gaussian with a FWHM of 15 Mm.
Figure 15 (to be compared with Fig. 11) demonstrates that the representation of the target function by the averaging kernel is indeed much better in the case of inversion using the Aztype target function and the p_{1} mode. A similar conclusion may be drawn when comparing Figs. 16 and 12 for the inversions with an ignored crosstalk. The superiority of the p_{1} inversion for δc_{s} is even more visible when comparing the horizontally averaged averaging kernels plotted as a function of depth (Fig. 17).
Fig. 15.
Averaging kernels for δc_{s} inversion at 0.5 Mm depth. The crosstalk was minimised. The Az target function was used. See Fig. 2 for details. 

Open with DEXTER 
Fig. 16.
Averaging kernels for δc_{s} inversion at 0.5 Mm depth. No crosstalk minimisation. The Az target function was used. See Fig. 2 for details. 

Open with DEXTER 
Fig. 17.
Left panel: cut through the line x = 0, y = 0 of the component of the averaging kernel and its comparison with the target function for the GdG fmode inversion. Right panel: same as in left panel but for the Az p_{1}mode inversion. 

Open with DEXTER 
8. Summary
The results shown in the figures presented here can also be represented by numbers. In Table 1 we give some statistical properties of the inversions for all physical quantities in question, that is, for three components of the flow velocity and the sound–speed perturbations.
For each quantity separately (by columns) we give the correlation coefficient of the inversion component in the direction of the inversion with the ideal (expected) answer, which is always positive and large. The correlation coefficient between the inverted estimate and the ideal answer is given, which is always positive, however, it is lower for the vertical flow and sound–speed perturbations. That is mainly due to the presence of the random noise. Finally, the RMS values of the component in the direction of the inversion, of the total crosstalk contribution, and of the random noise help to evaluate the S/N of the inverted estimates and the significance of the crosstalk pollution. The comparison between the upper and lower tables demonstrates the expected performance of the crosstalk minimisation.
In previous section, we performed an inversion for the sound–speed perturbations with a different target function and a different set of sensitivity kernels to improve the inversion. The improvements of such an approach are also apparent in Table 2 (that should be compared to the last column in Table 1). The correlation coefficient of the inverted estimate is much larger than before, the randomnoise level is smaller, and also the crosstalk pollution is smaller even in the case when it is not minimised.
Correlation coefficients and RMS of the inversions using the synthetic travel times.
Correlation coefficients and RMS of inversions using the synthetic travel times.
Moreover, in Fig. 18 we plot examples of the inversion weights for Δ ≈ 19 Mm. The weights are well localised around the central point.
Fig. 18.
Examples of the weight functions of the inversion with the minimisation of the crosstalk and the GdG target function. Rows correspond to individual inversions. Columns correspond to individual averaging geometries. One can see the contributions of the mean geometry to the vector flows, especially to the vertical component. 

Open with DEXTER 
9. Conclusions
We introduced the improved methodology of the SOLA time–distance inversions for local helioseismology. Our improvements consist of incorporating the mean traveltime geometry in the inversions and also of the use of the traveltime sensitivity kernels for all flow components and the sound–speed perturbations at once. This new methodology allows us to invert, for the first time, for the full value of the vertical velocity (not only for its variations as before) and to quantify the crosstalk between the vector flows and the sound–speed perturbations.
We validated our inversions by using the synthetic data based on the realistic numerical simulation of the solar convection provided by Matthias Rempel (Rempel 2014; DeGrave et al. 2014). The use of the synthetic data allowed us to properly quantify all the inversion components and to gain empirical experience of the use of the improved pipeline. To summarise the results:
For υ_{x} and υ_{y} (the horizontal flow) it is not necessary to minimise the crosstalk (at least near the surface). Our methodology essentially does not constitute an improvement to what has been used recently by other authors. We only point out that the targeted RMS of the random noise should be 20 m s^{−1} or less if one inverts the nearsurface flows using travel times averaged over a day or so.
– For υ_{z} (the vertical flow) it is essential to minimise the crosstalk, otherwise the results are nullified by the leakage mainly from the horizontal components. The leakage from the sound–speed perturbations is acceptably small. Our findings confirm the results seen in the past. The RMS of the random noise should be lower than 2 m s^{−1} for the nearsurface flows (again, for oneday averaged travel times). An obvious advantage of our new methodology is that one can also quantify the largescale (with scales larger than the typical horizontal extent of the averaging kernel) component of the vertical flow and not only its horizontal variations as before. The inversion is sensitive also to the largescale offset because of the incorporation of the mean traveltime kernels. Also, it seems that the incorporation of the mean travel times helps to regularise the crosstalk contribution. By using the improved methodology, we achieved significant improvements in inversions for the vertical flow.
At first sight it seems that the minimisation of the crosstalk in the case of the inversion for δc_{s} does not improve the results. This is because the leakage component from the other quantities is positively correlated with the sound–speed perturbations and thus the leakage “helps” to “measure” δc_{s}. Such an approach leads to an overestimation of the magnitude of the sound–speed variations by almost a factor of two in the nearsurface layers. Such an inversion, however, is not the inversion for δc_{s}. When the crosstalk minimisation is introduced, the agreement with the ideal answer decreases, which is in part due to the worse fit of the target function by the averaging kernel, and also by the lower S/N achieved for the travel times averaged over 24 h. The case we present here cannot be further improved because only the f mode, which is not particularly sensitive to the sound–speed perturbations, was used. The reliability improves by a large extent when the first acoustic ridge is used in a similar setup.
Acknowledgments
D.K. is supported by the Grant Agency of Charles University under grant No. 532217. M.Š. is supported by the project RVO:67985815 and the grant project 1806319S awarded by Czech Science Foundation. We are grateful to Matthias Rempel for providing us with the snapshot of the numerical simulation. The sensitivity kernels were computed by the code KC3 kindly provided by Aaron Birch. The authors would like to thank the anonymous referee for the valuable comments and suggestions, which greatly improved the quality of the paper. This research has made use of NASA Astrophysics Data System Bibliographic Services. Author contributions: M.Š. and D.K. designed the research. D.K. implemented and performed the research within his Ph.D projects under the supervision of M.Š. D.K. drafted the paper and both authors contributed to the final manuscript.
References
 Backus, G., & Gilbert, F. 1968, Geophys. J. Int., 16, 169 [NASA ADS] [CrossRef] [Google Scholar]
 Birch, A. C., & Gizon, L. 2007, Astron. Nachr., 328, 228 [NASA ADS] [CrossRef] [Google Scholar]
 Brüggen, M., & Spruit, H. C. 2000, Sol. Phys., 196, 29 [NASA ADS] [CrossRef] [Google Scholar]
 Brun, A. S., Miesch, M. S., & Toomre, J. 2004, ApJ, 614, 1073 [NASA ADS] [CrossRef] [Google Scholar]
 Burston, R., Gizon, L., & Birch, A. C. 2015, Space Sci. Rev., 196, 201 [NASA ADS] [CrossRef] [Google Scholar]
 ChristensenDalsgaard, J., Dappen, W., Ajukov, S. V., et al. 1996, Science, 272, 1286 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 DeGrave, K., Jackiewicz, J., & Rempel, M. 2014, ApJ, 788, 127 [NASA ADS] [CrossRef] [Google Scholar]
 Duvall, T. L., & Hanasoge, S. M. 2013, Sol. Phys., 287, 71 [NASA ADS] [CrossRef] [Google Scholar]
 Duvall, T. L., Hanasoge, S. M., & Chakraborty, S. 2014, Sol. Phys., 289, 3421 [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]
 Fournier, D., Gizon, L., Hohage, T., & Birch, A. C. 2014, A&A, 567, A137 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gizon, L., & Birch, A. C. 2004, ApJ, 614, 472 [NASA ADS] [CrossRef] [Google Scholar]
 Gizon, L., Birch, A. C., & Spruit, H. C. 2010, ARA&A, 48, 289 [NASA ADS] [CrossRef] [Google Scholar]
 Greer, B. J., Hindman, B. W., Featherstone, N. A., & Toomre, J. 2015, ApJ, 803, L17 [NASA ADS] [CrossRef] [Google Scholar]
 Hanasoge, S. M., Duvall, T. L., & Sreenivasan, K. R. 2012, Proc. Nat. Acad. Sci., 109, 11928 [NASA ADS] [CrossRef] [Google Scholar]
 Hansen, P. C. 1999, The Lcurve and its use in the Numerical Treatment of Inverse Problems (Lyngby, Denmark: IMM, Department of Mathematical Modelling, Technical University of Denmark) [Google Scholar]
 Jackiewicz, J., Birch, A. C., Gizon, L., et al. 2012, Sol. Phys., 276, 19 [NASA ADS] [CrossRef] [Google Scholar]
 Kosovichev, A. G., & Zhao, J. 2016, in Lecture Notes in Physics, eds. J.P. Rozelot, & C. Neiner (Berlin: Springer Verlag), 914, 25 [NASA ADS] [CrossRef] [Google Scholar]
 Louis, A. K., & Maass, P. 1990, Inverse Prob., 6, 427 [NASA ADS] [CrossRef] [Google Scholar]
 Miesch, M. S., Featherstone, N. A., Rempel, M., & Trampedach, R. 2012, ApJ, 757, 128 [NASA ADS] [CrossRef] [Google Scholar]
 Pijpers, F. P., & Thompson, M. J. 1992, A&A, 262, L33 [NASA ADS] [Google Scholar]
 Rempel, M. 2014, ApJ, 789, 132 [NASA ADS] [CrossRef] [Google Scholar]
 Rempel, M., Schüssler, M., & Knölker, M. 2009, ApJ, 691, 640 [NASA ADS] [CrossRef] [Google Scholar]
 Scherrer, P. H., Schou, J., Bush, R. I., et al. 2012, Sol. Phys., 275, 207 [NASA ADS] [CrossRef] [Google Scholar]
 Schou, J. 1997, in AAS/Solar Physics Division Meeting #28, Bull. Am. Astron. Soc., 29, 913 [Google Scholar]
 Schou, J., Scherrer, P. H., Bush, R. I., et al. 2012, Sol. Phys., 275, 229 [NASA ADS] [CrossRef] [Google Scholar]
 Tong, C. H., Thompson, M. J., Warner, M. R., & Pain, C. C. 2003, ApJ, 596, L263 [NASA ADS] [CrossRef] [Google Scholar]
 Švanda, M. 2012, ApJ, 759, L29 [NASA ADS] [CrossRef] [Google Scholar]
 Švanda, M. 2013, ApJ, 775, 7 [NASA ADS] [CrossRef] [Google Scholar]
 Švanda, M., Gizon, L., Hanasoge, S. M., & Ustyugov, S. D. 2011, A&A, 530, A148 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Švanda, M., Roudier, T., Rieutord, M., Burston, R., & Gizon, L. 2013, ApJ, 771, 32 [NASA ADS] [CrossRef] [Google Scholar]
 Zhao, J., & Kosovichev, A. G. 2004, ApJ, 603, 776 [NASA ADS] [CrossRef] [Google Scholar]
 Zhao, J., Georgobiani, D., Kosovichev, A. G., et al. 2007, ApJ, 659, 848 [NASA ADS] [CrossRef] [Google Scholar]
All Tables
Correlation coefficients and RMS of the inversions using the synthetic travel times.
Correlation coefficients and RMS of inversions using the synthetic travel times.
All Figures
Fig. 1.
Map of correlation between the vertical velocity and the sound–speed perturbation in the nearsurface layers obtained from the realistic simulation of the solar convection by Rempel (2014). This plot is analogous to Fig. 8 of Švanda et al. (2011). The black circle corresponds to the width of the target function (15 Mm). 

Open with DEXTER  
In the text 
Fig. 2.
Averaging kernels for υ_{x} inversion at 0.5 Mm depth. The crosstalk was minimised. The red curve corresponds to the halfmaximum of the GdG target function, the yellow curve corresponds to the halfmaximum of the averaging kernel and the blue solid and blue dotted (not present in this figure) curves correspond to +5% and −5% of the maximum of the averaging kernel, respectively. In the top row there are horizontal slices of the averaging kernel at the target depth and in the bottom row there are vertical slices perpendicular to the symmetries. In the right part one can see the selection of the tradeoff parameters and the standard deviation when assuming the travel times averaged over 24 hours. 

Open with DEXTER  
In the text 
Fig. 3.
Similar to Fig. 2: averaging kernels for υ_{x} inversion at 0.5 Mm depth without crosstalk minimisation. 

Open with DEXTER  
In the text 
Fig. 4.
Top row: inversion with minimisation of the crosstalk. Bottom row: inversion without minimisation of the crosstalk. First column: contribution to . Second column: contribution to . Third column: contribution to . Fourth column: contribution to . Fifth column: noise contribution to . Sixth column: . Seventh column: ideal answer. 

Open with DEXTER  
In the text 
Fig. 5.
Averaging kernels for υ_{z} inversion at 0.5 Mm depth. The crosstalk was minimised. See Fig. 2 for details. 

Open with DEXTER  
In the text 
Fig. 6.
Averaging kernels for υ_{z} inversion at 0.5 Mm depth. No crosstalk minimisation. See Fig. 2 for details. 

Open with DEXTER  
In the text 
Fig. 7.
Horizontal integral of the mean geometry contribution of various physics quantities to the averaging kernel. 

Open with DEXTER  
In the text 
Fig. 8.
Individual contributions into and its comparison with the ideal answer. See Fig. 4 for details. 

Open with DEXTER  
In the text 
Fig. 9.
Top row: with minimisation of the crosstalk and without minimisation of the crosstalk. Bottom row: and the ideal answer. The correlation coefficient between and the ideal answer is 0.87. The colour bars are not the same. 

Open with DEXTER  
In the text 
Fig. 10.
Top row: results of the inversion for υ_{z} with the mean and the difference geometries. Middle row: results of the inversion with the difference geometries only. Bottom row: difference of two aforementioned inversions. Left column: . Right column: ideal answers of the inversions. We point out that the colour bars are not scaled. 

Open with DEXTER  
In the text 
Fig. 11.
Averaging kernels for δc_{s} inversion at 0.5 Mm depth. The crosstalk was minimised. See Fig. 2 for details. 

Open with DEXTER  
In the text 
Fig. 12.
Averaging kernels for δc_{s} inversion at 0.5 Mm depth. No crosstalk minimisation. See Fig. 2 for details. 

Open with DEXTER  
In the text 
Fig. 13.
Individual contributions from and noise to and its comparison with the ideal answer. See Fig. 4 for details. 

Open with DEXTER  
In the text 
Fig. 14.
Top row: with minimisation of the crosstalk and without minimisation of the crosstalk. Bottom row: and ideal answer. The correlation coefficient between and ideal answer is 0.90. The colour bars are not the same. 

Open with DEXTER  
In the text 
Fig. 15.
Averaging kernels for δc_{s} inversion at 0.5 Mm depth. The crosstalk was minimised. The Az target function was used. See Fig. 2 for details. 

Open with DEXTER  
In the text 
Fig. 16.
Averaging kernels for δc_{s} inversion at 0.5 Mm depth. No crosstalk minimisation. The Az target function was used. See Fig. 2 for details. 

Open with DEXTER  
In the text 
Fig. 17.
Left panel: cut through the line x = 0, y = 0 of the component of the averaging kernel and its comparison with the target function for the GdG fmode inversion. Right panel: same as in left panel but for the Az p_{1}mode inversion. 

Open with DEXTER  
In the text 
Fig. 18.
Examples of the weight functions of the inversion with the minimisation of the crosstalk and the GdG target function. Rows correspond to individual inversions. Columns correspond to individual averaging geometries. One can see the contributions of the mean geometry to the vector flows, especially to the vertical component. 

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.