Open Access
Issue
A&A
Volume 634, February 2020
Article Number A44
Number of page(s) 16
Section The Sun
DOI https://doi.org/10.1051/0004-6361/201937007
Published online 04 February 2020

© B. Proxauf et al. 2020

Licence Creative Commons
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Open Access funding provided by Max Planck Society.

1. Introduction

Recently, Löptien et al. (2018, hereafter LGBS18) discovered global-scale Rossby waves in maps of flows on the surface of the Sun. These waves are waves of radial vorticity that may exist in any rotating fluid body. Even though Rossby waves were predicted to exist in stars more than 40 years ago (Papaloizou & Pringle 1978; Saio 1982), solar Rossby waves were difficult to detect because of their small amplitudes (∼1 m s−1) and long periods of several months. Solar Rossby waves contain almost as much vorticity as large-scale solar convection. The dispersion relation of solar Rossby waves is close to the standard relation for sectoral modes, ω = −2Ω/(m + 1), where Ω is the rotation rate of a rigidly rotating star and m is the azimuthal order (Saio 1982). Rossby waves have a retrograde phase speed and a prograde group speed. In LGBS18, the authors also measured the horizontal eigenfunctions, which peak at the equator.

The detection of solar Rossby waves was confirmed by Liang et al. (2019, hereafter LGBD19) with time-distance helioseismology (Duvall et al. 1993) using data covering more than 20 years, obtained from the Solar and Heliospheric Observatory (SOHO) and from the Solar Dynamics Observatory (SDO; Pesnell et al. 2012). Alshehhi et al. (2019), in an effort to speed up ring-diagram analysis (RDA; Hill 1988) via machine learning, also saw global-scale Rossby waves. Hanasoge & Mandal (2019) and Mandal & Hanasoge (2019) provide another recent Rossby wave confirmation using a different technique of helioseismology known as normal-mode coupling (Woodard 1989; Hanasoge et al. 2017).

Knowledge about the latitude dependence of Rossby wave eigenfunctions is incomplete, as LGBS18 studied only their real parts. In a differentially rotating star, the horizontal eigenfunctions are not necessarily spherical harmonics (and may not even separate in latitude and depth). Also, little is known observationally about the depth dependence of the Rossby waves. It would be well worth distinguishing between the few existing theoretical models of the depth dependence (Provost et al. 1981; Smeyers et al. 1981; Saio 1982; Wolff & Blizard 1986).

In this paper, we explore the latitude dependence of the eigenfunctions, as well as the phase and amplitude of solar Rossby waves as functions of depth from the surface down to more than 16 Mm using helioseismology. We use observations from the Helioseismic and Magnetic Imager (HMI; Schou et al. 2012) on board SDO, processed with RDA. From these we attempt to measure the eigenfunctions of the Rossby waves in the solar interior. For comparison near the surface, we also use data from local correlation tracking of granulation (LCT; November & Simon 1988).

2. Data and methods

We used maps of the horizontal velocity, derived from two different techniques applied to SDO/HMI observations. The first dataset consists of LCT (granulation tracking) flow maps at the surface (Löptien et al. 2017) and covers almost six years from May 20, 2010 to March 30, 2016. The second dataset comprises RDA flow maps from the HMI ring-diagram pipeline (Bogart et al. 2011a,b; see also Bogart et al. 2015). For comparisons with LCT, we took a period as close to the LCT period as possible, i.e., May 19, 2010 to March 31, 2016, while for all other results we used a longer period of more than seven years from May 19, 2010 to December 29, 2017; this corresponds to 102 Carrington rotations (CRs), i.e., CR 2097 - 2198.

2.1. Overview of LCT data

The LCT flow maps are obtained from and processed as described in Löptien et al. (2017). They are created by applying the Fourier LCT code (FLCT; Welsch et al. 2004; Fisher & Welsch 2008) to track the solar granulation in pairs of consecutive HMI intensity images. The image pairs are separated by 30 min. Several known systematic effects such as the shrinking-Sun effect (Lisle et al. 2004; Löptien et al. 2016) and effects related to the SDO orbit are present in the LCT maps. Therefore the maps are decomposed into Zernike polynomials, a basis of 2D orthogonal functions on the unit disk, and the time series of the coefficient amplitudes for the lowest few Zernike polynomials are filtered to remove frequencies of one day and one year (associated with the SDO orbit) as well as all harmonics up to the Nyquist frequency. The zero frequency is also removed. The filtered maps are then tracked at the sidereal Carrington rate and remapped onto an equi-spaced longitude-latitude grid with a step size of 0.4° in both directions.

2.2. Overview of ring-diagram data

The ring-diagram pipeline (Bogart et al. 2011a,b) takes HMI Dopplergrams as input and remaps them onto tiles spanning 182 × 182 Mm (i.e., 15° each in latitude λ and longitude φ at the equator). The tiles overlap each other by roughly 50% in each direction such that the tile borders fall onto the centers of adjacent tiles. Both the latitude and longitude sampling are half the tile size. The latitude grid is linear and includes the equator, while the longitude grid is also linear, but is latitude-dependent. Each tile is tracked for 1728 min (28.8 hr) at the sidereal Carrington rate. The temporal grid spacing is, on average, 1/24 of the synodic Carrington rotation period of 27.2753 days.

In the pipeline, for each tile a 3D local power spectrum is computed from the tracked Dopplergrams. The velocity fit parameters Ux,  nℓ (prograde) and Uy,  nℓ (northward) are extracted via a ring-fit algorithm (Haber et al. 2000) for different solar oscillation modes, which are indexed by their radial order n and angular degree . The flow velocities ux and uy are inferred for various target depths via a 1D optimally localized averages (OLA) inversion. The inversion results for the six-parameter fits of the 15° tiles sample a range of target depths from 0.97 R to 1 R (step size 0.001 R), corresponding to a nonlinear grid of measurement depths (median of the ring-diagram averaging kernels) from 0.976 R to 1 R. In this paper, the term depth always refers to measurement depth and not to target depth.

The inversion results are stored in the Joint Science Operations Center (JSOC) data series hmi.V_rdvflows_fd15_frame. However, up to inversion module rdvinv v.0.91, the inversion results depended on the input tile processing order due to an array initialization bug. This caused significantly lower velocity uncertainties for tiles near latitude 7.5° and Stonyhurst longitude 37.5°, even when averaged over seven years, but also slightly affected the velocities. At the same disk locations the bug caused a correlation of ux with the B0 angle. Since rdvinv v.0.92 is officially only applied since March 2018, we re-inverted the entire dataset ourselves for the analysis shown in this work.

Apart from the array initialization bug, we found several other issues with the default HMI ring-diagram pipeline that have not yet been solved. Among these are under-regularization in the inversion for some individual tiles, leading to relatively narrow averaging kernels and anomalously high noise. Finally, the number of ring fits used for the inversion depends strongly on disk position. This may lead to systematic effects and additional noise.

The ring-diagram velocities ux reported at a certain measurement depth r at the equator for an angular rotation rate Ω(r) are equal to Ω(r)R instead of the local velocity Ω(r)r. Since we are interested in the latter, we multiplied ux by r/R. By analogy, we also applied this factor to uy and to all other latitudes. Additionally, the inversion does not account for the quantity βnℓ, defined, for example, in Eq. (3.357) of Aerts et al. (2010). The quantity βnℓ is related to the effect of the Coriolis force on the mode frequency splitting. For uniform rotation in particular, at fixed m, βnℓ completely describes the effect of the rotation on the mode frequency splitting. Both issues are described in more detail in Appendix A.

2.3. Post-processing of ring-diagram data

The ring-diagram data are organized in CRs, which undergo several processing steps, including the removal of systematic effects, an interpolation in longitude, an interpolation in time, and the removal of limb data.

Several systematic effects are present in the ring-diagram velocities, such as center-to-limb effects that depend on the disk position of the tile (Baldner & Schou 2012; Zhao et al. 2012). There are time-independent effects and systematics with a one-year period, which are probably related to the B0 angle. To remove the systematics, we fit the time series at each position on the disk (in Stonyhurst coordinates) with sinusoids

(1)

and subtract the fits from the flow velocities. We used all available CRs to determine the fit parameters.

Because of the specific tile coordinate selection used by the ring-diagram pipeline (Bogart et al. 2011a), which seeks to optimally cover the visible disk, tile centers at different latitudes have Stonyhurst longitudes that are offset by multiples of 2.5° from each other. To obtain a latitude-independent longitude grid, we interpolated the flow maps in Stonyhurst longitude using splines (Appendix B).

We also interpolated the ring-diagram flows in time similarly with splines to fill missing time steps due to instrumental issues (only 12 out of 2448 time steps), which cause a too low observational duty cycle (≤​70%). We interpolated the data in the Carrington reference frame so as to use always roughly the same physical locations on the Sun. This mixes different systematics, which are primarily dependent on disk position, but we should already have removed the dominant contributions at this stage. We interpolated every missing time step from roughly the same number of data points (all available time steps within the corresponding disk passage) using splines (Appendix B).

The output uncertainties from the ring-diagram pipeline increase strongly toward the limb, in particular beyond an angular great-circle distance of roughly 65° to the crossing of the central meridian with the equator (λ = 0° ,φ = 0°). We thus only used ring-diagram data within 65° of (λ = 0° ,φ = 0°).

2.4. From velocity maps to power spectra of radial vorticity

From this stage onward ring-diagram and LCT data are processed similarly. The processing steps include a shift to the equatorial rotation rate νeq = Ωeq/2π = 453.1 nHz, the subtraction of the longitude mean, the calculation of the radial vorticity, a spherical harmonic transform (SHT), and a Fourier transform of the SHT coefficient time series.

The flow maps are shifted from the tracking rate (sidereal Carrington rate) to the surface sidereal equatorial rotation rate of 453.1 nHz, an average of zonal flows inferred from global-mode analysis of SDO/HMI observations (Larson & Schou 2018). We shifted the LCT data in Fourier space via a time-dependent phase factor, applying the same convention for the Fourier transform as LGBS18. The ring-diagram data are first apodized by a raised cosine in angular great-circle distance to the point (λ = 0° ,φ = 0°) to suppress near-limb data and are shifted via spline-interpolation (Appendix B).

We next subtracted the longitude mean from the data to remove any remaining large-scale flows. Differential rotation and meridional circulation should have already been subtracted in the RDA or LCT post-processing, but any possible longitude-independent flows still in the data are removed in this step.

Subsequently, we calculated the radial vorticity (via second-order central finite differences) as follows:

(2)

where r is the measurement depth. We decomposed the resulting maps into spherical harmonics and performed a temporal Fourier transform of the spherical harmonic coefficients. Last, we calculated the power and the phase (where the phase range is the half-open interval (−180° ,180°]). The sign convention is such that waves with positive m and negative frequency ν have a retrograde phase speed.

If not stated otherwise, the terms power spectrum or Fourier transform used in this paper always refer to the power spectrum or Fourier transform of the radial vorticity. Similarly, we discuss eigenfunctions of radial vorticity. These eigenfunctions are not spherical harmonics, however (LGBS18). In particular, while the modes can be meaningfully indexed by m, the angular degree is not observable. Throughout the paper thus only refers to the projection of the Rossby wave modes onto the corresponding spherical harmonic and not to the Rossby wave eigenfunction itself. We also use the terms latitudinal and radial eigenfunctions, which assumes separability in the r and λ coordinates. This assumption is addressed in more detail in Sect. 4.

3. Results

3.1. Radial vorticity maps

Figure 1 shows example vorticity maps from LCT surface flows and from RDA flows near the surface, at intermediate, and at large depths (0.7, 9.9, and 16.5 Mm), averaged over the first rotation in the dataset (May 20, 2010 to June 16, 2010). The LCT data have a much better horizontal resolution than the ring-diagram data and thus pick up small-scale convective contributions. To be able to compare LCT with RDA, we thus smooth the LCT vorticity with 1D Gaussian filters of width σ = 6° both in latitude and longitude.

thumbnail Fig. 1.

Radial vorticity maps from LCT at the surface and from RDA at depths 0.7, 9.9, and 16.5 Mm. The radial vorticity is averaged over one rotation (from May 20, 2010 to June 16, 2010). The LCT map is smoothed in latitude and longitude with a Gaussian filter (σ = 6°) to filter out small-scale convection.

Open with DEXTER

We do not expect perfect agreement of the two methods because of their different sensitivities to horizontal scales and to different depths. Nonetheless, the LCT map shows similar features as the near-surface (0.7 Mm) ring-diagram map. While large absolute radial vorticities are visible at high latitudes (beyond ±50°) in the LCT but not in the ring-diagram data, the vorticities near the equator agree. As a test, we interpolate the LCT data to the RDA grid using a 2D bicubic spline. The correlation coefficient between the interpolated LCT and the ring-diagram maps decreases with the latitude width of the strip of pixels considered and there is a steep decrease beyond ±45°. The correlation is 0.92 when including only equatorial pixels, 0.79 for pixels within ±45°, and 0.59 for all pixels, i.e., within ±52.5°. The noise increases strongly with depth (see lower panels of Fig. 1), but the main vorticity features are still visible.

3.2. Power spectra of radial vorticity

Figure 2 shows the Rossby wave power of the  = m component for the ring-diagram data near the surface (0.7 Mm) versus frequency and azimuthal order m (LGBS18 detected only the sectoral component of the Rossby waves). We divide the power, at each m, by the frequency average over [−300, 100] nHz near the surface (0.7 Mm). The visible power ridge corresponds to the Rossby wave signal. The mode frequency increases with m roughly according to the textbook dispersion relation for sectoral waves, ω = −2Ωeq/(m + 1), as seen earlier by LGBS18.

thumbnail Fig. 2.

Sectoral power spectrum ( = m) of the radial vorticity for RDA data at depth 0.7 Mm. The solid red line indicates the Rossby wave frequencies from LGBD19 for m = 3 and from LGBS18 for m ≥ 4. Frequency intervals for the Rossby wave region and the background region are indicated by the red and black dashed lines, respectively. The ridge of power at positive frequencies is due to the first side lobe of the window function. For better visibility of the Rossby waves at low m, the power is normalized at each m by the frequency average over [−300, 100] nHz. The color scale is truncated at 50% of the maximum value (black).

Open with DEXTER

Besides the Rossby wave signal there are other ridges, that, at fixed Δm = m − m′, are shifted from the Rossby waves by roughly Δm(νeq−1/(1 yr)), where νeq − 1/(1 yr)∼421.4 nHz. The main contribution to these side lobes comes from a temporal window function, which is introduced by solar rotation and not by time gaps; the time coverage is very good (see Sect. 2). This leads to side lobes of wave power from modes at m′ to modes at m. We only show the side lobes for Δm = +1, but we typically observe the side lobes above the noise between Δm = −2 and Δm = +3. In LGBD19, the authors use 21 years of time-distance data from a combined sample of observations from the Michelson Doppler Imager (MDI) on board SOHO and from SDO/HMI and they discuss the window function in detail.

Figure 3 shows the power versus frequency for different azimuthal orders m. We divide the power, at each m, by the frequency average of the  = m = 8 mode over [−300, 100] nHz near the surface (0.7 Mm). The power decreases from 0.7 to 9.9 Mm, then increases toward 16.5 Mm, but the depth dependence is modest (≤​​​20%). We also see that the wave power decreases with m faster for RDA than for LCT owing to the different sensitivity kernels, as found by LGBS18. The  = m = 3 signal has a multi-peak structure and is thus difficult to measure. We do not observe Rossby waves for  = m ≤ 2; the dash-dotted blue lines for  = m = 1 and  = m = 2 indicate the expected mode frequencies from the textbook dispersion relation.

thumbnail Fig. 3.

Sectoral power spectra ( = m) showing the Rossby wave power in the LCT data (black line) and the RDA data at depths 0.7 and 16.5 Mm (cyan and red lines). The power is normalized by the average of the m = 8 power in the range [−300, 100] nHz. The dash-dotted vertical lines in the m = 1 and m = 2 panels indicate the frequencies ω = −2Ωeq/(m + 1). The frequency axes of the m = 1 and m = 2 panels are shifted with respect to the other panels.

Open with DEXTER

The wave power side lobes due to the window function explain why the  = m = 6 side lobe in Fig. 2 even exceeds the main signal: the adjacent  = m = 7 mode has a higher relative power (see Fig. 3). Systematic effects that are fixed in the Stonyhurst reference frame can be easily misinterpreted as an  = m = 1 Rossby wave signal (see the LCT curve in Fig. 3), as their frequency (the rotation rate) is equal to the  = m = 1 Rossby wave frequency.

We assume that there is background power contributing to the observed power at the Rossby peak, but measuring its contribution directly at the peak is impossible. Since we are limited by the side lobes, we use a region halfway between the peak and the next side lobe, i.e., shifted from the peak by half the rotation rate. We checked that the shift direction does not matter much, so for the central background frequency, we use the Rossby wave frequencies from LGBD19 and LGBS18 for m = 3 and m ≥ 4, respectively, plus half the rotation frequency νeq. We use the full widths at half maximum γref = Γref/2π from LGBD19 and LGBS18 for m = 3 and m ≥ 4, respectively, and perform a least-squares second-order polynomial fit in m to obtain a smoothed linewidth γsmooth. We use γsmooth for the width of the peak and background frequency intervals. Thus our peak and background frequency intervals at each m are

(3)

These definitions are used in the analysis of latitudinal eigenfunctions in Sect. 3.3.

In Fig. 2, we see that the peak interval (dashed red lines) typically captures the main wave power well. The 1D power spectra, however, reveal that the background interval (not shown in Fig. 3, but see dashed black lines in Fig. 2), however, is potentially contaminated by scattered signal power, for example, for  = m = 5, 6, and 14. To check how the frequency interval definition affects our results, we performed our analysis for several different peak and background intervals. The results are consistent, thus we adopt Eq. (3) for the peak and background intervals.

Unlike LGBS18, we see evidence for non-sectoral components of the Rossby waves. For  = m + 2 the 2D power spectrum shows, for 6 ≤ m ≤ 13, a ridge of power at very similar frequencies to those of the  = m Rossby waves seen in Fig. 2, apart from a higher relative noise level and side lobes. This is confirmed by the 1D cuts at fixed values of m. We do not see structure in the power spectra for  = m + k other than for k = 0 and k = 2. In Sect. 3.3.3, we indeed show that the latitudinal eigenfunctions of Rossby waves are not sectoral spherical harmonic functions (in agreement with LGBS18).

3.3. Latitudinal eigenfunctions of Rossby waves

To estimate the latitudinal eigenfunctions, we first remove small-scale convection from the LCT maps via smoothing with a 6° Gaussian in latitude. Next we compute the Fourier transform of the radial vorticity maps ζ(t, r, λ, φ) from LCT and RDA in time and longitude as follows:

(4)

The variables are discrete and take values at time steps tj = jT/Nt (integer 0 ≤ j <  Nt), longitudes φk = 2πk/Nφ (integer 0 ≤ k <  Nφ), frequencies νs = s/T (integer s, with −Nt/2 ≤ s ≤ Nt/2 − 1 for even Nt), and azimuthal orders m (integer, with −Nφ/2 ≤ m ≤ Nφ/2 − 1 for even Nφ). In this case, T, Nt, and Nφ are the observation period and the number of data points in time and longitude, respectively. We apply a filter to select the Rossby waves one m at a time, i.e.,

(5)

The filter Fm(ν) is equal to one within the Rossby wave ridge and zero elsewhere. Since ζ(t, r, λ, φ) is real, the symmetry applies. We then transform back to time to obtain

(6)

In this way we obtain filtered time-latitude vorticity maps for every m. Because there is no symmetry , the filtered vorticity maps are in general complex.

LGBS18 do a similar analysis for LCT data, in particular for rotation-averaged maps and filtering within ±30 nHz around the central mode frequencies. We do the entire latitudinal eigenfunction analysis for LCT and RDA, for full time-resolution maps and maps averaged in time within individual solar rotations, and for a ±27 nHz (five frequency pixels) and a linewidth filter (Eq. (3)) around the central mode frequencies. The different time-resolution and filtering cases yield consistent results; we thus show only the outcome for the full time-resolution and linewidth filtering. However, LGBS18 take the real part of the complex . This is equivalent to assuming that the phase of the eigenfunction is independent of latitude. We address the implications of this in Sect. 3.3.3 in more detail.

To estimate uncertainties for all results in this paper, we split the data into equal-size time intervals, apply our analysis to each chunk, and calculate the standard deviation over the results (for complex quantities separately for the real and imaginary part). Appendix C gives more details on error estimation and validation. Because of the small number of chunks, low-number statistics are an issue and the reported error bars are relatively uncertain.

For the sake of clarity, for the simple case of a single m Rossby wave with a single frequency νm and an eigenfunction Cm(r, λ), the vorticity field would be given by

(7)

We apply two different methods to obtain the eigenfunctions Cm(r, λ) near the surface, the covariance method (Sect. 3.3.1), and the SVD method (Sect. 3.3.2). The former is used also by LGBS18.

3.3.1. Covariance

We calculate, at each m, the temporal covariance of the vorticity between the equator near the surface (target depth R = R − 0.7 Mm for RDA) and all other latitudes and depths, normalized by the variance at the equator near the surface

(8)

where the angle brackets ⟨ ⋅ ⟩t denote a temporal average and is the centered vorticity. The function Cm(r, λ) is complex-valued since is in general complex. By construction Cm(r = R, λ = 0° ) is unity. Appendix D shows that Cm can also be obtained by a linear fit to the vorticity. The same covariance can be computed with the LCT data.

3.3.2. Singular value decomposition

We present a second, new method to obtain latitudinal eigenfunctions. We want to separate the filtered vorticity at each azimuthal order m and depth r, i.e., a 2D matrix, into a latitude and a time dependence, i.e.,

(9)

Applying a singular value decomposition (SVD), we can decompose the vorticity as

(10)

where s(r,  m), j is the singular value of index j with left and right singular vectors U(r,  m), j and V(r,  m), j and k is the minimum between the number of grid points in time and latitude. The square of s(r,  m), j measures the variance captured by its singular vectors. By convention the singular values are sorted in descending order, thus the first singular vector contains more variance than any other individual singular vector.

Assuming that there is only one nonzero singular value, s(r,  m), 0, the SVD gives the desired decomposition of the vorticity into one time and one latitude function. Our observations indeed have one clearly dominant singular value.

Given that the noise at high latitudes increases steeply, we crop our vorticity maps for the SVD to latitudes within ±50° of the equator. Also, the SVD does not account for the varying noise of the remaining latitudes. To ensure that latitudes with larger uncertainties are given less weight, we filter the original vorticity maps once more in Fourier space for the noise, calculate the temporal standard deviation σm of the noise-filtered maps, and compute . We filter for the noise by taking either all frequencies except for five pixels around the peak or all frequencies within the background interval (see Eq. (3)). The two different filters give consistent results. At each m, the SVD is performed on the weighted maps and the resulting latitude vectors are multiplied by σm again to undo the weighting. We apply the weighting only to LCT, since the ring-diagram data are already apodized (see Sect. 2). We select the first latitude singular vector near the surface and normalize it by its value at the equator.

3.3.3. Results for the latitudinal eigenfunctions

Figures 4 and 5 show the real and imaginary parts of the horizontal eigenfunctions of Rossby waves versus latitude for different m. The real part is consistent with the findings from LGBS18. The imaginary part, however, was not discussed by LGBS18.

thumbnail Fig. 4.

Real part of Cm(λ) for different azimuthal orders m and four different methods (see legend in panel m = 4). The shaded areas indicate the 1σ error estimates.

Open with DEXTER

thumbnail Fig. 5.

Imaginary part of Cm(λ) for different azimuthal orders m and four different methods (see legend in panel m = 4). The shaded areas indicate the 1σ error estimates. For comparison, the solid gray curves show the real part of Cm for the LCT covariance-based data. The plotting ranges are the same as in Fig. 4.

Open with DEXTER

In the current paper, we find that the LCT and the RDA results are mostly consistent for the near-surface layers. Also, almost all m show agreement between the covariance and SVD results. This in particular holds for the modes with the largest amplitudes, i.e., for 7 ≤ m ≤ 10. On the other hand, the modes m = 4 and to a lesser extent m = 15, where Rossby wave measurements become difficult, display larger errors but nonetheless consistent results. The m = 3 results for the different techniques disagree and are noisy. The m = 5 and m = 6 results for the real part differ slightly between the covariance and SVD methods. While the covariance yields a real part of the eigenfunction quite similar to those of other modes, the SVD-based results show maxima around latitudes of ±10 − 15°. Apparently, there the SVD picks up some variance that is uncorrelated with the equator. It is unclear whether it is just noise, or a real signal of a different kind of latitudinal eigenfunctions.

The eigenfunction shape is similar for different modes. The real part decreases away from the equator, flips sign, and approaches zero after going through a local minimum. The imaginary part is much noisier than the real part, as indicated by the error estimates. For most m, it is close to zero and flat near the equator, but reaches minima at high latitudes. The latitude of the minima appears to move toward the equator with increasing m.

As can be seen from, for example, the red curves in Fig. 5, the imaginary part appears to be mostly positive for 3 ≤ m ≤ 6. For 7 ≤ m ≤ 9 the sign of the imaginary part is unclear. For 10 ≤ m ≤ 15, the imaginary part is predominantly negative. The presence of an imaginary part induces a phase for the latitudinal eigenfunctions that can be interpreted as a latitude-dependent shift of the sinusoid in longitude. A positive sign of the imaginary part means that the horizontal eigenfunctions at high latitudes are leading in the retrograde direction with respect to the equator. Conversely, a negative sign would indicate that the eigenfunctions at high latitudes are trailing with respect to the equator. This may provide important constraints on the theory of latitudinal eigenfunctions of Rossby waves.

Figure 4 suggests that the real part of the eigenfunctions is more confined to low latitudes for higher values of m. We study the m-dependence of several characteristic parameters illustrated in Fig. 6, i.e., the width at an eigenfunction real part of Re(C) = 0.5 (a half width at half maximum; HWHM), the latitude of the eigenfunction real part sign reversal, λ0, and the latitude and value of the minimum, λmin and Cmin, respectively. To reduce the noise level we derive the eigenfunctions from maps symmetrized in latitude before measuring these parameters.

thumbnail Fig. 6.

Schematic description of the real part of Cm(λ) for a given m. The various parameters that describe the curve are the HWHM, the latitude at zero crossing (λ0), the latitude at minimum (λmin), and the minimum value (Cmin).

Open with DEXTER

To obtain the latitude widths at Re(C) = 0.5 and Re(C) = 0, we linearly fit the two points closest to these Re(C) values. The latitude and value of the minimum are obtained by quadratically fitting three points around the minimum derived without fitting. We do this to avoid oscillating RDA results due to the coarse 7.5° latitude sampling. For LCT, the effects of fitting the minimum (or not) are minimal. There are no results for m = 3 and m = 4 because of the poor quality and different shape of the eigenfunctions. We already stated the difficulties in characterizing these modes. As described at the beginning of Sect. 3.3 and in Appendix C, to derive uncertainties, we compute the standard deviation over the results for different time chunks, separately for the real and the imaginary part.

Table 1 shows how these parameters, measured for the LCT data from the covariance method, depend on m. Although not given in the table, we also measure the parameters for the RDA and the SVD results. We thus also discuss the m-dependence for those measurements; this dependence is mostly consistent with that of the LCT covariance-based parameters.

Table 1.

Parameters of the real part of Cm(λ) for the LCT covariance-based data; see Fig. 6.

The latitude width at Re(C) = 0.5 indeed decreases with m, quasi-linearly between m = 7 and m = 13. The slope is roughly −1° per m. The decrease might flatten off at high m, but this could also be caused by noise. We observe slightly different latitude widths between the covariance and SVD eigenfunctions at low m for Re(C) = 0.5, but similar widths at Re(C) = 0. Toward higher m, λ0 is consistent with a flat profile, until around m = 13 the eigenfunction widths become smaller. The latitude of the minimum, λmin, shows an m-(in)dependence similar to λ0. There is a strong discrepancy for m = 13 between LCT and RDA, indicating that this mode is not trivial to characterize. This could be caused by noise. To some extent we could already see this in the power spectrum in Fig. 3, where the m = 13 linewidth is large compared to all other m. The error on λmin might be underestimated here, since as seen in the asymmetric eigenfunctions in Fig. 4 the minimum of the LCT data is more poorly defined for m = 13 than for other modes. Finally, the value of the minimum, Cmin, is different between the different analysis methods at m = 5 and m = 6, as seen before. Otherwise, it is quasi m-independent and has at most a very mild increase with m from m = 7 onward, which is likely covered by noise, however.

As mentioned before, the latitudinal eigenfunctions appear to have two nodes (zero crossings) at latitudes ±λ0. This in combination with Fig. 2 and the subsequent discussion indicates that the eigenfunctions have significant contributions from  = m and  = m + 2 components. To quantify these contributions, we project the symmetric eigenfunctions Cm(λ) onto associated Legendre polynomials , to obtain the coefficients

(11)

The sum goes over all latitudes λ = /Nλ (integer −Nλ/2 ≤ k <  Nλ/2), where Nλ is the number of data points in latitude. The are normalized such that . The associated Legendre polynomials used in the decomposition are not orthogonal over the limited observed latitude range. However, we do not expect this to be a problem since we see later in this section that the near-sectoral associated Legendre polynomials, whose amplitude is concentrated toward the equator, are the dominant contributions to the latitudinal eigenfunctions. Because of the symmetry of the eigenfunctions, only cℓm with even  − m ≥ 0 are nonzero. We find that the real and the imaginary parts of the eigenfunctions for almost all m can be approximated well (within 1σ) when using only the contributions from cℓm for m ≤  ≤ m + 6, except for m = 3, which is very noisy. The approximation also does not work well at the high latitudes (beyond ±40°) for the real part (for some modes) and at the near-equatorial latitudes for the imaginary part (for almost all modes).

Table 2 shows the coefficients cℓm for m ≤  ≤ m + 6 for the LCT covariance-based latitudinal eigenfunctions. As usual the uncertainties are calculated from the standard deviation over the coefficients for different time chunks (Appendix C), separately for the real and the imaginary part. The real part of the eigenfunctions is clearly dominated by the  = m component. The contribution from the  = m + 2 component is significant as well and has a relative strength of 30−50%. This is consistent with the observations from the 2D and 1D power spectra in Sect. 3.2. The real part of the cmm and cm+2, m each depend weakly on m. The real part of several of the coefficients with larger is insignificant. The imaginary part, on the other hand, has significant, dominant contributions at  = m + 4 for m ≤ 10 and at  = m + 2 for m ≥ 11, whereas the  = m and  = m + 6 components are often insignificant. The term insignificant refers to an absolute value of cℓm of less than 1σ. Nonetheless, independent of the estimated error bars, 11 out of 12 modes within 4 ≤ m ≤ 15 have the same sign for cm + 4,  m, suggesting that the  = m + 4 contribution to the real part is significant, despite being within 1σ from zero. A similar argument holds for the imaginary part of the latitude dependence of the Rossby wave eigenfunctions.

Table 2.

Coefficients cℓm for the LCT covariance-based data.

For the latitudinal eigenfunctions of Rossby waves there are so far only a few theoretical studies such as Lee & Saio (1997) and Townsend (2003). These studies typically gave either analytic (asymptotic) expressions and/or numeric calculations, but the former expressions do not agree well with the latter calculations for Rossby waves (Townsend 2003). Although both studies indicate that the latitudinal eigenfunctions are not concentrated near the equator, we cannot sensibly compare their findings to our measurements. In particular these models assume a uniform rotation rate. Also these authors used the traditional approximation, i.e., they neglected the horizontal component of the rotation vector. This approximation requires the squared Brunt-Väisälä frequency N2 to be much higher than both the squared oscillation frequency ω2 and the squared rotation rate Ω2. The validity of the traditional approximation thus has to be critically examined within the convection zone of the Sun.

3.4. Radial eigenfunctions of Rossby waves

3.4.1. Depth-dependent ring-diagram systematics

To study the Rossby wave depth dependence, we must check to which depths RDA is reliable. For this we compare the solar rotation profile from ring-diagram velocities with the results from SDO/HMI global modes from the JSOC data series hmi.V_sht_2drls (Larson & Schou 2018). The global modes have a 72-day time sampling from April 30, 2010 to June 4, 2017, a 1.875° latitude sampling, and a nonlinear depth grid with many more points near the surface than at larger depths. Global modes are expected to give a precise and accurate solar rotation profile. We interpolate the global mode results to the ring-diagram latitude-depth grid via 2D bicubic splines, which is reasonable because the global-mode inversions do not vary on scales of their original grid; we then average the 72-day chunks over time. The chunk scatter of the rotation rate is used to estimate the uncertainty. We divide the ring-diagram data into five intervals of length 480 time steps (20 rotations), average the chunks over time and estimate the error from the scatter, convert the velocities into rotation rates, and add the sidereal Carrington rate to correct for the ring-diagram tracking.

Figure 7 shows the equatorial rotation rate versus depth from global modes and ring-diagram velocities, both averaged over longitude and at Stonyhurst longitudes of ±52.5° (the outermost longitudes in our vorticity maps). The global modes yield a smooth profile with extremely small errors. The ring-diagram data show a small offset at small depths, but, more importantly, inconsistency with the global modes at large depths. Of course, it is difficult to judge how well the results should agree because of the different kernels of the datasets and thus different depth (and latitude) sensitivities. The −52.5° longitude curve has a small local maximum around 8 Mm. For longitudes even further east (not shown), the rotation rate has an even stronger excess (a bump) there.

thumbnail Fig. 7.

Solar equatorial rotation rate as a function of depth. The global-mode helioseismology result is given by the black curve. The blue curve is for RDA after averaging over all longitudes. The green and red curves show the ring-diagram results at Stonyhurst longitudes ±52.5°. The shaded areas give the 1σ error estimates. The observations cover more than seven years. The dashed black lines indicate the depth range within which the ring-diagram results are in best agreement with each other.

Open with DEXTER

The most worrisome point is the disagreement between different ring-diagram longitudes themselves and also with the longitude average, below roughly 8 Mm (indicated by the left dashed line in Fig. 7). Because we averaged the data over more than seven years, any short-lived flows and even longer-lived structures should be filtered away and the longitude gradient from east to west should thus not exist. This points to a deeper problem with the ring fits and the pipeline processing that generated these fits. The presence of systematic effects in HMI ring-diagram data has also been extensively discussed in Komm et al. (2015).

Finally, we note that Fig. 7 is affected by an issue related to the ring-diagram inversion, since the inversion does not account for the quantity βnℓ. A discussion of this issue and a brief check of the magnitude of the effect is given in Appendix A. The latter showed that the main effect is a depth-independent underestimation of the ring-diagram velocities by 1 − 2 m s−1 or equivalently of the rotation rate by less than 0.5 nHz. This does not affect our main conclusions. The small, but significant difference between the rotation rates from global modes and ring diagrams cannot be caused by the βnℓ issue (it has the wrong sign), but may possibly instead be due to different averaging kernel widths, systematics, or other unknown effects.

3.4.2. Determining the Rossby wave depth dependence

In this section, we discuss only the sectoral ( = m) component of the power spectrum of radial vorticity. The Rossby wave power Pm(ν, r) and phase Φm(ν, r) thus depend on frequency, depth, and azimuthal order. Based on the assumption of damped oscillations and stochastic wave excitation, we perform a maximum-likelihood Lorentzian fit (Anderson et al. 1990) to the power spectra for the longer ring-diagram period, separately at each m. We use the functional form

(12)

We fit all the depths (except for the surface, i.e., r = 0.0 Mm, where the ring-diagram data are unreliable) at once, with a common central frequency ν0,  m and linewidth γm, but with individual amplitudes Am(r) and backgrounds Bm(r). The Lorentzian fit of the power spectra, in most cases, fits well to the observations. As seen in Fig. 3, the  = m = 6 and  = m = 13 modes have large linewidths and their power spectra show fine structure. The  = m = 3 mode has been fit by LGBD19, but not by LGBS18.

To determine error bars for the amplitudes and backgrounds via chunked data (Appendix C), we also perform the Lorentzian fit for each chunk separately, fitting again all depths together, but keeping the central frequency and linewidth fixed at the fit results of ν0,  m and γm from Eq. (12) for the full time period (to prevent unstable fits). Because we keep these parameters fixed for each chunk, we cannot derive their uncertainties based on the standard deviation over the chunks. We thus do a Monte Carlo simulation and generate 1000 realizations of synthetic power spectra according to Eq. (C.1) (Appendix C) and perform the Lorentzian fit for each realization analogously to the fit for the observations. The median of the parameters over the Monte Carlo realizations is consistent with the fit parameters for the observations. While the error based on the Monte Carlo simulation contains realization noise, the model we use (Lorentzian and stochastic excitation) does not include all features of the observed power spectra. The chunk-based error likely describes the physical system more accurately, by also including other variance contributions, such as from temporal effects on the Rossby waves; we could imagine, for example, solar cycle effects. This may also explain the discrepancy between the two types of errors of order 30% (Appendix C). This disagreement is, however, small enough to not affect the significance of the results for the radial eigenfunctions. We thus use the uncertainties based on the Monte Carlo simulation for the central frequency and the linewidth.

Table 3 compares the fit parameters from this study with the results from LGBD19 and LGBS18. As in LGBD19 and LGBS18, the upper and lower errors give the difference between the quantiles comprising the central 68.3% (the distributions are non-Gaussian) and the fit parameters for the observations. We also calculate the uncertainties on the central frequency following Libbrecht (1992) and find that those (symmetric) errors typically underestimate the Monte Carlo quantile errors by roughly 1 nHz. A possible reason for this could be that we use a finite frequency fitting interval.

Table 3.

Measured frequencies and linewidths of the Rossby waves from RDA sectoral power spectra with azimuthal orders in the range 3 ≤ m ≤ 15.

The fit parameters for the observations and those from LGBD19 and LGBS18 typically agree within 1σ or better. The central frequencies and the linewidths for the  = m = 6 mode differ by 10 and 37 nHz, but the fit is sensitive to the fitting range. The  = m = 3 fit parameters do not agree. In LGBD19, the authors, using 21 years of data, observed that the multi-peak structure of the  = m = 3 power spectrum (Fig. 3) seen in data with shorter periods collapses to a narrow single peak, which indicates that the discrepancy of the fit parameters is explained by stochastic excitation of the Rossby waves and vanishes when fitting data with a longer time period. Our errors are typically more symmetric and often smaller than those of LGBD19 and LGBS18. The lower errors for the linewidth often agree better than the upper errors, indicating a tail of high values (skewness) in the LGBD19 and LGBS18 estimate distributions. Reasons for the differences in the error estimates may lie in the simultaneous fitting of all depths at once or in the different observation periods of our datasets and those of LGBS18.

To determine the power depth dependence, we use the amplitude A derived from the Lorentzian fit (see Eq. (12)) and we define the normalized power of the signal as

(13)

We thus normalize by the depth average of the amplitude of the Lorentzian. The normalized power is independent of temporal amplitude variations due to Rossby wave excitation.

3.4.3. Results for the radial eigenfunctions

Figure 8 shows the depth dependence of the  = m = 8 phase, but the behavior is similar for other m. For easier comparison, we remove phase jumps of 360° and move the depth average to zero. The phase at the frequency of maximum power is almost constant with depth, within roughly ±3°. The phase at the background, at the center of the background interval, varies much more strongly with depth, within roughly ±100°, although the phases at other background frequencies sometimes show much less variation. The background phase is not random in depth (see Appendix C). In particular, phase changes are gradual and smooth; the depths are correlated. This could indicate a significant contribution from scattered signal power to the background. Nonetheless, peak and background display distinctly different depth dependences. We also find that phases at different frequencies across the peak and background are different. Frequency averages of phases are thus not useful. However, as seen, for single frequencies the phase at the peak is nearly constant with depth, while the background phase varies with depth.

thumbnail Fig. 8.

Phase at a frequency of roughly −87.4 nHz corresponding to the peak of power for the Rossby mode  = m = 8, as a function of depth (blue curve). The green line shows the phase of the background at the center of the background interval (see Eq. (3)). The depth refers to the median of the ring-diagram averaging kernels (orange dots), which corresponds to certain target depths (black dots).

Open with DEXTER

Figure 8 also shows the main parameters of the ring-diagram averaging kernels for a few target depths, i.e., the first, second (median) and third quartiles and the width (interquartile range). The flow measurements are well-localized near the surface, but smeared out over a broad depth range at large depths. The ring-diagram depth covariance matrix (not shown) indicates a similar behavior and shows that different depths are mostly independent near the surface, while at the largest depths there is high correlation and they thus do not give independent results. This could maybe also explain why the background phase is not random in depth. At large depths, the center of the averaging kernels (second quartile) moves away from the target depth, but the averaging kernels are relatively symmetric.

The background power for different m (not shown) generally increases with depth and at least for some modes there could be a minimum at 8 − 9 Mm, albeit with little significance given the large errors.

Figure 9 shows the signal power (Eq. (13)). The quantity 𝒫signal typically decreases from the surface toward a depth of 8 Mm, significantly as shown by the errors. Even further inside the Sun the power often increases again and reaches near-surface or even higher values. The 1σ errors shown in this plot give the standard deviation, but they do not indicate 68.3% probability intervals, since the power distribution is non-Gaussian (power cannot be negative). More information about error estimation can be found in Appendix C.

thumbnail Fig. 9.

Blue lines show the Rossby wave power 𝒫signal as a function of depth for different values of m. The blue shaded areas indicate the 1σ errors. The orange curves are fits of the form const. × r2α over depths between 0.7 and 7.4 Mm (between the vertical dashed lines). The orange dashed curves are extrapolations to larger depths.

Open with DEXTER

Provost et al. (1981) presented a theoretical argument that the Rossby wave eigenfunctions for the horizontal displacement are proportional to rm under the assumption that the modes are incompressible. Thus, under this theory, the radial vorticity is expected to be proportional to rm−1. To compare this to our observations, we perform a fit of the form const. × r2α to 𝒫signal within the dashed black lines (0.7 to 7.4 Mm) where the RDA is more reliable (see Fig. 7). We assume that the data points are uncorrelated in depth. Obviously, the fit does not reproduce the increase of power at large depths.

Figure 10 compares the observed and theoretical exponent α. The fitted exponent has very large error bars. It is consistent with the theoretical model from Provost et al. (1981), but also with the absence of any trend with m. Although the exponent depends strongly on the fit range because of the kink at roughly 5 Mm in Fig. 9, we also do not find inconsistency with a flat dependence on m within other fit intervals. Thus the current error estimates do not allow a definitive statement on the radial dependence of Rossby waves.

thumbnail Fig. 10.

Exponent α as a function of m, measured in the top 7.4 Mm (blue line) and 1σ error (blue shaded area). The dashed black line corresponds to the model α = m − 1, obtained under the assumption of non-divergent motions (Provost et al. 1981).

Open with DEXTER

4. Summary

We build on LGBS18, who investigated Rossby waves mostly using granulation tracking, by studying several Rossby wave properties via the analysis of radial vorticities computed from RDA at different depths and LCT at the surface. We obtained several new results: independently the latitudinal eigenfunctions with RDA (including a more complete, complex-valued description of the eigenfunctions), and the Rossby wave power and phase depth dependence.

We calculated latitudinal eigenfunctions of Rossby waves from the radial vorticity maps via the covariance between the equator and different latitudes and from the singular vectors of an SVD. We confirmed the shape of the real part of the eigenfunction from LGBS18, who used the covariance method on symmetrized LCT data. We also saw consistency between covariance and SVD results, except for m = 5 and m = 6, where the SVD eigenfunctions had maxima around ±10−15° instead of at the equator. The shape of the real part of the latitudinal eigenfunctions seen for most m indicates that the Rossby waves have maximum amplitudes near the equator, as found by LGBS18. The imaginary part appears to be mostly positive for low m (3 ≤ m ≤ 6); this part varies around zero for intermediate m (7 ≤ m ≤ 9) and is mostly negative for high m (10 ≤ m ≤ 15). A nonzero imaginary part may be due to attenuation and to the interaction of the waves with large-scale flows. In particular, the interaction of viscous Rossby waves with latitudinal differential rotation leads to the formation of critical layers at intermediate latitudes (Gizon & Fournier, priv. comm.).

We defined and measured characteristic parameters for the real part of the eigenfunctions and we found that the width at an eigenfunction value of 0.5 (the HWHM) decreased with m, in contrast to the m-independent width at a value of 0 and the latitude and value of the eigenfunction minimum. We also decomposed the eigenfunctions into associated Legendre polynomials and saw that the real part is dominated by  = m and  = m + 2 contributions, while the imaginary part consists mostly of  = m + 4 and  = m + 2 contributions for low and high m, respectively.

We compared rotation rates from ring-diagram data and global modes as functions of depth and saw a small offset at small depths and disagreement at large depths, but most importantly inconsistency of different ring-diagram longitudes. This indicated systematic effects in the ring-diagram pipeline (see also Komm et al. 2015).

We studied the Rossby wave power and phase depth dependence in detail for the first time. The phase at the peak is stable with depth, in contrast to the phase of the background. The background power almost monotonically increases with depth, while the signal power decreases toward a depth of 8−9 Mm and then increases again. The radial eigenfunctions of the Rossby waves are (at small depths) consistent with a power-law decrease, in particular both with the theoretical Provost et al. (1981) model (exponent m − 1) and an m-independent exponent. However, the Provost et al. (1981) model is based on assumptions that are not exactly correct for the Sun (e.g., uniform rotation). We can constrain the radial dependence of the eigenfunctions only very weakly owing to the high uncertainties on the observed exponents.

The analysis presented in this paper implicitly makes the assumption that the Rossby wave eigenfunctions are separable in depth and latitude. Our data show a similar latitude dependence for the different depths, separability thus appears to be a reasonable assumption. The results shown in this work motivate further research on Rossby wave eigenfunctions, which is a necessary condition for the interpretation of the measured mode frequencies.

Acknowledgments

B. Proxauf is a member of the International Max Planck Research School for Solar System Science at the University of Göttingen; he conducted the data analysis and contributed to the interpretation of the results and to the writing of the manuscript. We thank Z.-C. Liang for providing help with the fitting of the modes in frequency space and V. Böning for providing beta values for Appendix A. The HMI data are courtesy of NASA/SDO and the HMI Science Team. The data were processed at the German Data Center for SDO funded by the German Aerospace Center DLR. We acknowledge partial support from the European Research Council Synergy Grant WHOLE SUN #810218.

References

  1. Aerts, C., Christensen-Dalsgaard, J., & Kurtz, D. W. 2010, Asteroseismology (Springer) [Google Scholar]
  2. Alshehhi, R., Hanson, C. S., Gizon, L., & Hanasoge, S. 2019, A&A, 622, A124 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  3. Anderson, E. R., Duvall, Jr., T. L., & Jefferies, S. M. 1990, ApJ, 364, 699 [NASA ADS] [CrossRef] [Google Scholar]
  4. Baldner, C. S., & Schou, J. 2012, ApJ, 760, L1 [NASA ADS] [CrossRef] [Google Scholar]
  5. Bogart, R. S., Baldner, C., Basu, S., Haber, D. A., & Rabello-Soares, M. C. 2011a, J. Phys. Conf. Ser., 271, 012008 [Google Scholar]
  6. Bogart, R. S., Baldner, C., Basu, S., Haber, D. A., & Rabello-Soares, M. C. 2011b, J. Phys. Conf. Ser., 271, 012009 [Google Scholar]
  7. Bogart, R. S., Baldner, C. S., & Basu, S. 2015, ApJ, 807, 125 [NASA ADS] [CrossRef] [Google Scholar]
  8. Christensen-Dalsgaard, J. 2008, Ap&SS, 316, 113 [NASA ADS] [CrossRef] [Google Scholar]
  9. Christensen-Dalsgaard, J. 2011, Astrophysics Source Code Library [record ascl:1109.002] [Google Scholar]
  10. Duvall, Jr., T. L., Jefferies, S. M., Harvey, J. W., & Pomerantz, M. A. 1993, Nature, 362, 430 [NASA ADS] [CrossRef] [Google Scholar]
  11. Fisher, G. H., & Welsch, B. T. 2008, in Subsurface and Atmospheric Influences on Solar Activity, eds. R. Howe, R. W. Komm, K. S. Balasubramaniam, & G. J. D. Petrie, ASP Conf. Ser., 383, 373 [NASA ADS] [Google Scholar]
  12. Haber, D. A., Hindman, B. W., Toomre, J., et al. 2000, Sol. Phys., 192, 335 [NASA ADS] [CrossRef] [Google Scholar]
  13. Hanasoge, S., & Mandal, K. 2019, ApJ, 871, L32 [NASA ADS] [CrossRef] [Google Scholar]
  14. Hanasoge, S. M., Woodard, M., Antia, H. M., Gizon, L., & Sreenivasan, K. R. 2017, MNRAS, 470, 1404 [NASA ADS] [CrossRef] [Google Scholar]
  15. Hill, F. 1988, ApJ, 333, 996 [NASA ADS] [CrossRef] [Google Scholar]
  16. Komm, R., González Hernández, I., Howe, R., & Hill, F. 2015, Sol. Phys., 290, 1081 [NASA ADS] [CrossRef] [Google Scholar]
  17. Larson, T. P., & Schou, J. 2018, Sol. Phys., 293, 29 [NASA ADS] [CrossRef] [Google Scholar]
  18. Lee, U., & Saio, H. 1997, ApJ, 491, 839 [NASA ADS] [CrossRef] [Google Scholar]
  19. Liang, Z.-C., Gizon, L., Birch, A. C., & Duvall, T. L. 2019, A&A, 626, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Libbrecht, K. G. 1992, ApJ, 387, 712 [NASA ADS] [CrossRef] [Google Scholar]
  21. Lisle, J., & Toomre, J. 2004, in SOHO 14 Helio- and Asteroseismology: Towards a Golden Future, ed. D. Danesy, ESA Special Publication, 559, 556 [NASA ADS] [Google Scholar]
  22. Löptien, B., Birch, A. C., Duvall, T. L., Gizon, L., & Schou, J. 2016, A&A, 590, A130 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Löptien, B., Birch, A. C., Duvall, T. L., et al. 2017, A&A, 606, A28 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Löptien, B., Gizon, L., Birch, A. C., et al. 2018, Nat. Astron., 2, 568 [NASA ADS] [CrossRef] [Google Scholar]
  25. Mandal, K., & Hanasoge, S. 2019, ArXiv e-prints [arXiv:1908.05890] [Google Scholar]
  26. November, L. J., & Simon, G. W. 1988, ApJ, 333, 427 [NASA ADS] [CrossRef] [Google Scholar]
  27. Papaloizou, J., & Pringle, J. E. 1978, MNRAS, 182, 423 [NASA ADS] [CrossRef] [Google Scholar]
  28. Pesnell, W. D., Thompson, B. J., & Chamberlin, P. C. 2012, Sol. Phys., 275, 3 [NASA ADS] [CrossRef] [Google Scholar]
  29. Provost, J., Berthomieu, G., & Rocca, A. 1981, A&A, 94, 126 [NASA ADS] [Google Scholar]
  30. Saio, H. 1982, ApJ, 256, 717 [NASA ADS] [CrossRef] [Google Scholar]
  31. Schou, J., Scherrer, P. H., Bush, R. I., et al. 2012, Sol. Phys., 275, 229 [NASA ADS] [CrossRef] [Google Scholar]
  32. Smeyers, P., Craeynest, D., & Martens, L. 1981, Ap&SS, 78, 483 [NASA ADS] [CrossRef] [Google Scholar]
  33. Townsend, R. H. D. 2003, MNRAS, 340, 1020 [NASA ADS] [CrossRef] [Google Scholar]
  34. Welsch, B. T., Fisher, G. H., Abbett, W. P., & Regnier, S. 2004, ApJ, 610, 1148 [NASA ADS] [CrossRef] [Google Scholar]
  35. Wolff, C. L., & Blizard, J. B. 1986, Sol. Phys., 105, 1 [NASA ADS] [CrossRef] [Google Scholar]
  36. Woodard, M. F. 1989, ApJ, 347, 1176 [NASA ADS] [CrossRef] [Google Scholar]
  37. Zhao, J., Nagashima, K., Bogart, R. S., Kosovichev, A. G., & Duvall, Jr., T. L. 2012, ApJ, 749, L5 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Issues of the ring-diagram inversions

In this Appendix, we discuss two issues regarding the ring-diagram pipeline inversions. In order to obtain local velocities at a certain measurement depth r, the reported pipeline velocities must be multiplied by r/R. Additionally the pipeline inversion does not take the quantity βnℓ (see, e.g., Aerts et al. 2010, Eq. (3.357)) into account and thus the reported inversion velocities ux are slightly incorrect.

To see this, we study a simple case. For now, let us assume that the ring diagrams are not tracked. The frequency perturbation δωnℓm of the mode indexed by radial order n, angular degree , and azimuthal order m due to a radial differential rotation rate Ω(r) is

(A.1)

where Knℓ is the normalized rotation kernel for that mode, i.e., (see, e.g., Aerts et al. 2010, Eq. (3.358)). On the other hand, ring diagrams assume that the velocity mode fits Ux,  nℓ are equal to a radial integral over the true velocity flow field ux(r) weighted by flow sensitivity kernels. Based on inspection of the pipeline, we think that the used HMI kernels are normalized rotation kernels Knℓ from Eq. (A.1). Thus

(A.2)

To connect the two equations in a simple case, consider the Doppler shift of a sectoral ( = m) mode as seen by a ring diagram at the equator, i.e.,

(A.3)

In this equation, kx is the wavenumber in the prograde direction, which is related to m via kx = m/R. We conclude that

(A.4)

This is not consistent with Eq. (A.2) since βnℓ is missing from Eq. (A.2). Additionally we see that ux(r) should be interpreted as RΩ(r) and not as the local linear velocity rΩ(r).

To see what happens if the tracking rate is not zero, we now suppose that we track at rotation rate ΩT. Equation (A.4) then becomes

(A.5)

where is the ring measurement in the rotating frame. We now define the local deviation from the tracking rate δΩ(r) = Ω(r)−ΩT. Then we obtain

(A.6)

The first term is the same form as in Eq. (A.4), while the second term is an offset that depends on n and .

The conversion factor r/R is multiplied onto the data before any analysis is performed for this paper; see Sect. 2.2. The offset due to βnℓ depends on the set of mode ring fits, but it should be mostly time-independent, since the ring-diagram mode set does not vary much with time. Thus the time-dependent Rossby waves should not be sensitive to this effect and the only affected result in this paper should be the comparison of rotation rates in Fig. 7.

To estimate the effect of βnℓ on the inversion result, for a given input flow ux we generate artificial ring fits via Eq. (A.6), on which we run the ring-diagram inversion module to retrieve the output velocities. To compute , we assume a depth-independent flow equivalent to the tracking rate (sidereal Carrington rate), thus δΩ(r) = 0. We thus check only the second term of Eq. (A.6) and neglect that βnℓ also appears in the first term as a scaling factor. However, the effect due to the second term should be much larger than that due to the first term, as ΩT is much larger than δΩ(r) for the ring diagrams.

We use βnℓ values provided by V. Böning (priv. comm.). These were computed from eigenfunctions obtained from the Aarhus adiabatic oscillation package (ADIPLS; Christensen-Dalsgaard 2008, 2011). We lose roughly 25% of the original ring-fit modes, as we only have βnℓ values up to frequencies of 5 mHz. However, this does not critically change the mode set used during the inversion. We replace the actual pipeline ring fits with the artificial data. We leave all other data, including uncertainties on mode-fit velocities, as is and perform the inversion. The aforementioned conversion factor of r/R is multiplied onto the output velocities ux.

We see that the effect of βnℓ does not depend much on depth and that the retrieved ux are on the order of only 1.5 m s−1 (equivalent to roughly 0.1% of the tracking rate). The reason for this is that the inversion gives much more weight to the high modes for which the uncertainties are comparatively small. These modes typically have βnℓ values around 0.999, thus 1 − βnℓ ∼ 0.1%. We performed this check exemplarily for a ring-diagram tile at the first time step in our dataset (May, 20, 2010) at the point (λ = 0° ,φ = 0°). However, tests using different tiles show that this result does not depend much on time or disk position.

The effect of the pipeline inversion not accounting for βnℓ is thus an underestimation of the true velocity fields by roughly 1 − 2 m s−1, or equivalently approximately 0.4 nHz. This difference would be visible in Fig. 7, but does not change our main conclusions.

Appendix B: Interpolation and apodization of ring-diagram velocities

We interpolate the ring-diagram velocities separately in time and longitude (see Sect. 2) with different functions, depending on the number of available data points:

(B.1)

Before we interpolate the ring-diagram velocities to the surface equatorial rotation rate, we apodize these velocities with a raised cosine in angular great-circle distance ρ to the point (λ = 0° ,φ = 0°), see Sect. 2), as follows:

(B.2)

where β defines the steepness of the raised cosine flanks. We choose β = 0.3. The quantity T defines the central position of the flanks. We choose T such that zero is reached at ρ = 67.5° (where there are no more valid pixels). Apodizing the ring-diagram velocities (with different β), or not, gives consistent results.

Appendix C: Error estimation and error validation

C.1. Error estimation via chunked data

The uncertainties on all results are derived by dividing the time series of vorticity maps (in total 2448 time steps, i.e., 102 rotations, for RDA) into equal-size chunks and calculating the scatter over the results for the chunks. We find that for chunks longer than a few months (roughly six rotations), the Rossby wave signature is visible in the power spectra. We make a compromise between noise level and chunk statistics and divide the dataset into five chunks that are 480 time steps long each (20 rotations).

For the latitudinal eigenfunctions, where we only study the shorter LCT period (i.e., 78 rotations), we first used four chunks of length 19 rotations (rotation-averaged maps) and 470 time steps (full time resolution, for RDA), but we obtained very large errors for the SVD method for different single m, where single chunks gave singular vectors different from the usual eigenfunction shape. We think that the noise in our filtered maps might have been detected as the dominant term in the decomposition. We thus use a chunking with three chunks of length 26 rotations and 625 time steps, where all chunks have the expected first singular vectors.

C.2. Error validation via Monte Carlo simulation

We validate the chunking approach for the depth dependence via a Monte Carlo simulation. As a plausible physical model for the Rossby wave power spectrum, we assume a Lorentzian profile and a background (constant in frequency), each with a χ2-distributed random variable (stochastic excitation). In analogy to Eq. (12) we generate 1000 realizations of synthetic data for the Fourier transform of the radial vorticity ℱsyn (not the power P), at each m, as

(C.1)

We fix the amplitude Am(r), background Bm(r), central frequency ν0,  m, and full width at half maximum γm via the fit parameters for the observations. Furthermore we assume that the random variable 𝒩A,  m(ν) is constant with depth, i.e., the signal is fully correlated in depth, while for the background, we take a random variable 𝒩B,  m(ν, r) that is uncorrelated in depth. The random variables are complex Gaussian variables (with independent real and imaginary parts) with zero mean and unit variance, independent for each frequency. We analyze the realizations in the same way as the observations.

We observe inconsistency between two Monte Carlo estimates (by roughly 30%): one from chunking (like the observations), averaged over the realizations, and the other from the scatter of the power over the realizations. We find that the discrepancy is due to the temporal correlation of the different chunks; the temporal correlation matrix has values around 0.1 on the first off-diagonal. The two quantities agree when using a weighted average with weights based on the temporal covariance matrix.

Both Monte Carlo estimates disagree with the error for the observed data. This could be due to the depth correlation of the observations. We determine a (noisy) estimate of the depth covariance and correlation matrix of the observed background power from the different background frequencies and find strong correlations between different depths, even between the largest depth and the surface, with values above 0.25 in the off-diagonal corners. We use the observed depth covariance matrix for the Fourier transform (not the power) to correlate the Gaussian background random variable of our realizations, via a Cholesky decomposition.

Correlating the background in depth noticeably improves the agreement between Monte Carlo and observed errors. However, there is still a remaining discrepancy that can be attributed to our model, which is missing some features of the observed power spectrum. We do not account for the window function nor a background decreasing with frequency present in the observations. A detailed analysis goes beyond the scope of this paper.

Appendix D: Relation of covariance to linear fit

The covariance method (Sect. 3.3.1) is conceptually equivalent to a linear fit of the vorticity at each depth r and latitude λ, i.e.,

(D.1)

Let us for simplicity assume that the vorticity is real. The latitude and depth dependence in this vorticity separation ansatz is contained in the fit parameter am(r, λ). Let us assume that only is uncertain. For zero-mean quantities (such as our vorticity maps ), the slope of a linear fit without intercept, i.e., am(r, λ), is given as

(D.2)

Assuming the time dependence is given by the surface equatorial vorticity time series, i.e., , in Eq. (8) we can identify am(r, λ) with Cm(r, λ). Equation (D.1) implies that am(r = R, λ = 0° ) is unity.

The main disadvantage of the covariance method is the assumption of a noise-free vorticity at the equator, , required so that the time-dependence fm(t) is noise-free and the vorticity is the only uncertain quantity of the fit.

All Tables

Table 1.

Parameters of the real part of Cm(λ) for the LCT covariance-based data; see Fig. 6.

Table 2.

Coefficients cℓm for the LCT covariance-based data.

Table 3.

Measured frequencies and linewidths of the Rossby waves from RDA sectoral power spectra with azimuthal orders in the range 3 ≤ m ≤ 15.

All Figures

thumbnail Fig. 1.

Radial vorticity maps from LCT at the surface and from RDA at depths 0.7, 9.9, and 16.5 Mm. The radial vorticity is averaged over one rotation (from May 20, 2010 to June 16, 2010). The LCT map is smoothed in latitude and longitude with a Gaussian filter (σ = 6°) to filter out small-scale convection.

Open with DEXTER
In the text
thumbnail Fig. 2.

Sectoral power spectrum ( = m) of the radial vorticity for RDA data at depth 0.7 Mm. The solid red line indicates the Rossby wave frequencies from LGBD19 for m = 3 and from LGBS18 for m ≥ 4. Frequency intervals for the Rossby wave region and the background region are indicated by the red and black dashed lines, respectively. The ridge of power at positive frequencies is due to the first side lobe of the window function. For better visibility of the Rossby waves at low m, the power is normalized at each m by the frequency average over [−300, 100] nHz. The color scale is truncated at 50% of the maximum value (black).

Open with DEXTER
In the text
thumbnail Fig. 3.

Sectoral power spectra ( = m) showing the Rossby wave power in the LCT data (black line) and the RDA data at depths 0.7 and 16.5 Mm (cyan and red lines). The power is normalized by the average of the m = 8 power in the range [−300, 100] nHz. The dash-dotted vertical lines in the m = 1 and m = 2 panels indicate the frequencies ω = −2Ωeq/(m + 1). The frequency axes of the m = 1 and m = 2 panels are shifted with respect to the other panels.

Open with DEXTER
In the text
thumbnail Fig. 4.

Real part of Cm(λ) for different azimuthal orders m and four different methods (see legend in panel m = 4). The shaded areas indicate the 1σ error estimates.

Open with DEXTER
In the text
thumbnail Fig. 5.

Imaginary part of Cm(λ) for different azimuthal orders m and four different methods (see legend in panel m = 4). The shaded areas indicate the 1σ error estimates. For comparison, the solid gray curves show the real part of Cm for the LCT covariance-based data. The plotting ranges are the same as in Fig. 4.

Open with DEXTER
In the text
thumbnail Fig. 6.

Schematic description of the real part of Cm(λ) for a given m. The various parameters that describe the curve are the HWHM, the latitude at zero crossing (λ0), the latitude at minimum (λmin), and the minimum value (Cmin).

Open with DEXTER
In the text
thumbnail Fig. 7.

Solar equatorial rotation rate as a function of depth. The global-mode helioseismology result is given by the black curve. The blue curve is for RDA after averaging over all longitudes. The green and red curves show the ring-diagram results at Stonyhurst longitudes ±52.5°. The shaded areas give the 1σ error estimates. The observations cover more than seven years. The dashed black lines indicate the depth range within which the ring-diagram results are in best agreement with each other.

Open with DEXTER
In the text
thumbnail Fig. 8.

Phase at a frequency of roughly −87.4 nHz corresponding to the peak of power for the Rossby mode  = m = 8, as a function of depth (blue curve). The green line shows the phase of the background at the center of the background interval (see Eq. (3)). The depth refers to the median of the ring-diagram averaging kernels (orange dots), which corresponds to certain target depths (black dots).

Open with DEXTER
In the text
thumbnail Fig. 9.

Blue lines show the Rossby wave power 𝒫signal as a function of depth for different values of m. The blue shaded areas indicate the 1σ errors. The orange curves are fits of the form const. × r2α over depths between 0.7 and 7.4 Mm (between the vertical dashed lines). The orange dashed curves are extrapolations to larger depths.

Open with DEXTER
In the text
thumbnail Fig. 10.

Exponent α as a function of m, measured in the top 7.4 Mm (blue line) and 1σ error (blue shaded area). The dashed black line corresponds to the model α = m − 1, obtained under the assumption of non-divergent motions (Provost et al. 1981).

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text 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 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.