Properties of solar Rossby waves from normal mode coupling and characterizing its systematics

TL;DR

This study analyzes 20 years of helioseismic data to characterize solar Rossby waves, revealing their properties and systematics, and introduces a novel method to address measurement challenges caused by spatial leakage.

Contribution

A new approach is developed to mitigate spatial leakage effects in measuring solar Rossby waves from helioseismic data.

Findings

Rossby-mode velocities are around 0.5 m/s at the surface.

The method successfully separates harmonic degrees despite spatial leakage.

Analysis of 20 years of data provides detailed Rossby wave properties.

Abstract

Rossby waves play an important role in mediating the angular momentum of rotating spherical fluids, creating weather on Earth and tuning exoplanet orbits in distant stellar systems (Ogilvie 2014). Their recent discovery in the solar convection zone provides an exciting opportunity to appreciate the detailed astrophysics of Rossby waves. Large-scale Rossby waves create subtle drifts in acoustic oscillations in the convection zone, which we measure using helioseismology to image properties of Rossby waves in the interior. We analyze 20 years of space-based observations, from 1999 to 2018, to measure Rossby-mode frequencies, line-widths, and amplitudes. Spatial leakage affects the measurements of normal model coupling and complicates the analysis of separating out specific harmonic degree and the azimuthal number of features on the Sun. Here we demonstrate a novel approach to overcome this difficulty and test it by performing synthetic tests. We find that the root-mean-square velocity of the modes is of the order of 0.5 m/s at the surface.

Tables

Table 1. Table 1: Analysis of HMI data. Mode frequency, ω 0 subscript 𝜔 0 \omega_{0} , amplitude A 𝐴 \sqrt{A} , full width at half maximum Γ Γ \Gamma , background power B 𝐵 B that give best Lorentzian fits to the observed B 𝐵 B -coefficient spectra are tabulated in the co-rotating frame. The second column in the table gives the theoretical frequencies of modes in a co-rotating frame with tracking frequency 453 453 453 nHz. For comparison, we list the observed frequencies from other two studies, Löptien et al. ( 2018 ) and Liang et al. ( 2018 ) . The fitted spectrum is plotted in Figure ( 6 ). Here, MC stands for mode coupling.

		MC	Löptien et al. (2018)	Liang et al. (2018)
s	$\frac{2\Omega }{s+1}$	${\omega }_0/\left(2\pi \right)$	${\omega }_0/\left(2\pi \right)$	${\omega }_0/\left(2\pi \right)$	$\Gamma /\left(2\pi \right)$	$\sqrt{A}$	B
	nHz	nHz	nHz	nHz	nHz	cm s-1
1	453	451$\pm$0.1	—	—	5$\pm$0.5	538$\pm$157	31$\pm$647
3	226.5	233$\pm$3	259	254$\pm$2	58$\pm$11	70$\pm$24	6$\pm$155
5	151	156$\pm$0.4	157$\pm$4	156$\pm$2	12$\pm$1	68$\pm$15	131$\pm$22
7	113	111 $\pm$0.1	112$\pm$4	110$\pm$4	6$\pm$1	124$\pm$66	119$\pm$19
9	90.6	76$\pm$4	86$\pm$6	${82}_{-5}^{+4}$	53$\pm$12	32$\pm$ 12	47$\pm$34
11	75.5	54$\pm$2	75$\pm$7	46$\pm$7	66$\pm$7	26$\pm$7	31$\pm$12
15	56.6	18$\pm$1	${47}_{-6}^{+7}$	${22}_{-3}^{+2}$	10$\pm$3	25$\pm$12	42$\pm$6

Table 2. Table 2: Same as in Table ( 1 ), except parameter values are obtained from analyzing MDI data. The fitted spectrum is shown in Figure ( 7 ).

		MC	Löptien et al. (2018)	Liang et al. (2018)
s	$\frac{2\Omega }{s+1}$	${\omega }_0/(2\pi$)	${\omega }_0/\left(2\pi \right)$	${\omega }_0/\left(2\pi \right)$	$\Gamma /\left(2\pi \right)$	$\sqrt{A}$	B
	nHz	nHz	nHz	nHz	nHz	cm s-1
1	453	451$\pm$ 0.1	—	—	5$\pm$1	391$\pm$ 136	1357$\pm$ 480
3	226.5	249$\pm$ 0.4	259	254$\pm$2	12$\pm$ 2	107 $\pm$ 37	306$\pm$ 72
5	151	153$\pm$ 1	157$\pm$4	156$\pm$2	15$\pm$ 3	38$\pm$ 14	97$\pm$ 22
7	113	112$\pm$ 0.6	112$\pm$4	110$\pm$4	14$\pm$ 2	55$\pm$ 15	93$\pm$ 24
9	90.6	85$\pm$ 3	86$\pm$6	${82}_{-5}^{+4}$	52$\pm$ 9	26$\pm$ 9	23$\pm$ 18
11	75.5	56$\pm$ 2	75$\pm$7	46$\pm$7	41$\pm$ 5	19$\pm$ 5	30$\pm$ 6
13	64.7	48$\pm$ 4	40$\pm$10	24$\pm$7	65$\pm$ 12	12$\pm$ 4	24$\pm$ 4
15	56.6	25$\pm$ 2	${47}_{-6}^{+7}$	${22}_{-3}^{+2}$	41$\pm$ 6	15$\pm$ 3	24 $\pm$ 4

Equations

$$B_{st}^{\sigma }(n,\ell )=\frac{{\sum }_{m,\omega }{H}_{\ell mst}^{\sigma }(\omega ){\varphi }_{\ell m}^{\ast }(\omega ){\varphi }_{\ell m+t}(\omega +\sigma +t\Omega )}{{\sum }_{m,\omega }|H_{\ell mst}^{\sigma }{|}^2},$$

$$H_{\ell mst}^{\sigma }(\omega )=-2\omega (-1{)}^{m+t}\sqrt{2s+1}\left(\begin{array}{ccc}\ell & s& \ell \\ -(m+t)& t& m\end{array}\right)L_{\ell m}^{\ell m}{L}_{\ell m+t}^{\ell m+t}{N}_{\ell }(R_{\ell m}^{\omega \ast }|R_{\ell m+t}^{\omega +\sigma +t\Omega }{|}^2+|R_{\ell m}^{\omega }{|}^2{R}_{\ell m+t}^{\omega +\sigma +t\Omega }),$$

$$R_{\ell m}^{\omega }=\frac{1}{({\omega }_{n\ell m}-i{\Gamma }_{n\ell }/2{)}^2-{\omega }^2},$$

$$u^{\sigma }(r,\theta ,\varphi )=\sum _{st}{w}_{st}^{\sigma }(r)\hat{\mathbf{r}}\times {\nabla }_h{Y}_{st}(\theta ,\varphi ),$$

$$B_{st}^{\sigma }(n,\ell )=N_{\ell st}^{\sigma }\sum _{{\ell }^{\prime },{\ell }^{\prime \prime },m,m^{\prime },s^{\prime },t^{\prime },\omega }{L}_{\ell m}^{{\ell }^{\prime },m^{\prime }}{L}_{\ell m+t}^{{\ell }^{\prime \prime }{m}^{\prime }+t^{\prime }}{\gamma }_{tm}^{\ell s\ell }{H}_{\ell \ell mt}^{\sigma \ast }{\gamma }_{t^{\prime }{m}^{\prime }}^{{\ell }^{\prime \prime }{s}^{\prime }{\ell }^{\prime }}{H}_{{\ell }^{\prime }{\ell }^{\prime \prime }{m}^{\prime }{t}^{\prime }}^{\sigma }{b}_{s^{\prime }{t}^{\prime }}^{\sigma }({\ell }^{\prime },{\ell }^{\prime \prime }),$$

$$N_{\ell st}^{\sigma }=\frac{1}{{\sum }_{m,\omega }|{\gamma }_{tm}^{\ell s\ell }{H}_{\ell \ell mt}^{\sigma }{|}^2},$$

$$b_{st}^{\sigma }(\ell ,{\ell }^{\prime })=f_{{\ell }^{\prime }-\ell ,s}{\int }_{\odot }dr{w}_{st}^{\sigma }(r){\mathcal{K}}_{n\ell }(r).$$

$$B_{s,-s}^{\sigma }(n,\ell )=N_{\ell st}^{\sigma }\sum _{{\ell }^{\prime },{\ell }^{\prime \prime },m,m^{\prime },s^{\prime }}{L}_{\ell m}^{{\ell }^{\prime },m^{\prime }}{L}_{\ell m-s}^{{\ell }^{\prime \prime }{m}^{\prime }-s^{\prime }}{\gamma }_{-sm}^{\ell s\ell }{H}_{\ell \ell m-s}^{\sigma \ast }{\gamma }_{-s{m}^{\prime }}^{{\ell }^{\prime \prime }{s}^{\prime }{\ell }^{\prime }}{H}_{{\ell }^{\prime }{\ell }^{\prime \prime }{m}^{\prime }-s^{\prime }}^{\sigma }{b}_{s^{\prime },-s^{\prime }}^{\sigma }({\ell }^{\prime },{\ell }^{\prime \prime }).$$

$$w_{st}^{\sigma }(r_0)=\sum _{n\ell }{\alpha }_{n\ell ;r_0}{B}_{st}^{\sigma }(n,\ell ),$$

$$\chi ={\int }_{\odot }dr12(r-r_0{)}^2{\left(\sum _{n,\ell }{\alpha }_{n\ell ;r_0}{f}_{0,s}{K}_{n\ell }(r)\right)}^2+\lambda \sum _{n\ell }{\mathcal{N}}_{n\ell }{\alpha }_{n\ell ;r_0}^2.$$

$$w_{st}^{\sigma }(r)=\sum _k{\beta }_{st}^{\sigma }{B}_k(r),$$

$$\chi =\sum _{n,\ell }{\frac{\left(B_{st}^{\sigma }(n,\ell )-f_{0,s}{\int }_{\odot }dr{w}_{st}^{\sigma }(r){\mathcal{K}}_{n\ell }(r)\right)}{{\mathcal{N}}_{n\ell }^{\sigma }}}^2+\lambda \sum {\left(\frac{d^2{w}_{st}^{\sigma }}{d{r}^2}\right)}^2,$$

$$\frac{1}{{\mathcal{N}}_{n\ell }^{\sigma }}{f}_{0,s}{\int }_{\odot }dr{w}_{st}^{\sigma }(r){\mathcal{K}}_{n\ell }(r)=\frac{1}{{\mathcal{N}}_{n\ell }^{\sigma }}{B}_{st}^{\sigma }(n,\ell ),$$

$$\lambda \frac{d^2{w}_{st}^{\sigma }(r_k)}{d{r}^2}=0,$$

$$v_s(\omega )=\sum _{s^{\prime }}{L}_s^{s^{\prime }}{w}_{s^{\prime }}\delta (\omega -({\sigma }_{t^{\prime }}+t^{\prime }\Omega )).$$

$${\sigma }_t=\frac{2\Omega }{2|t|+1},$$

$$B_{s,-s}^{\sigma }(n,\ell )=\int dr\sum _{s^{\prime }}{\Theta }_s^{s^{\prime }}(n,\ell ,\sigma ,r)w_{s^{\prime }-s^{\prime }}^{\sigma }(r),$$

$${\Theta }_s^{s^{\prime }}(\sigma )=N_{\ell s-s}^{\sigma }\mathcal{K}(n,\ell )(r)\sum _{{\ell }^{\prime },{\ell }^{\prime \prime },m,m^{\prime },\omega }{L}_{\ell m}^{{\ell }^{\prime },m^{\prime }}{L}_{\ell m-s}^{{\ell }^{\prime \prime }{m}^{\prime }-s^{\prime }}{\gamma }_{-s^{\prime }m}^{\ell s\ell }{H}_{\ell \ell m-s}^{\sigma \ast }{\gamma }_{-s^{\prime }{m}^{\prime }}^{{\ell }^{\prime \prime }{s}^{\prime }{\ell }^{\prime }}{H}_{{\ell }^{\prime }{\ell }^{\prime \prime }{m}^{\prime }-s^{\prime }}^{\sigma },$$

$$B_{s,-s}^{\sigma }(n,\ell )=\int dr{\Theta }_s^s(n,\ell ,\sigma ,r)w_{s,-s}^{\sigma }(r).$$

$$F(\omega )=\frac{A}{1+[(\omega -{\omega }_0)/(\tau /2){]}^2}+B,$$

$$B_{st}^{\sigma }(n,\ell )=N_{\ell st}^{\sigma }\sum _{{\ell }^{\prime },{\ell }^{\prime \prime },m,m^{\prime },s^{\prime },t^{\prime },\omega }{L}_{\ell m}^{{\ell }^{\prime },m^{\prime }}{L}_{\ell m+t}^{{\ell }^{\prime \prime }{m}^{\prime }+t^{\prime }}{\gamma }_{tm}^{\ell s\ell }{H}_{\ell \ell mt}^{\sigma \ast }{\gamma }_{t^{\prime }{m}^{\prime }}^{{\ell }^{\prime \prime }{s}^{\prime }{\ell }^{\prime }}{H}_{{\ell }^{\prime }{\ell }^{\prime \prime }{m}^{\prime }{t}^{\prime }}^{\sigma }{b}_{s^{\prime }{t}^{\prime }}^{\sigma }({\ell }^{\prime },{\ell }^{\prime \prime }).$$

$$B_{st}^{\sigma }(n,\ell )=N_{\ell st}^{\sigma }\sum _{{\ell }^{\prime },{\ell }^{\prime \prime },m,m^{\prime },s^{\prime },t^{\prime },\omega }{\delta }_{\ell {\ell }^{\prime }}{\delta }_{m{m}^{\prime }}{\delta }_{\ell {\ell }^{\prime \prime }}{\delta }_{m+t{m}^{\prime }+t^{\prime }}{\gamma }_{tm}^{\ell s\ell }{H}_{\ell \ell mt}^{\sigma }{\gamma }_{t^{\prime }{m}^{\prime }}^{{\ell }^{\prime \prime }{s}^{\prime }{\ell }^{\prime }}{H}_{{\ell }^{\prime }{\ell }^{\prime \prime }{m}^{\prime }{t}^{\prime }}^{\sigma }{b}_{s^{\prime }{t}^{\prime }}^{\sigma }({\ell }^{\prime },{\ell }^{\prime \prime })$$

$$=N_{\ell st}^{\sigma }\sum _{m{s}^{\prime }\omega }{\gamma }_{tm}^{\ell s\ell }{\gamma }_{tm}^{\ell s^{\prime }\ell }|H_{\ell \ell mt}^{\sigma \ast }{|}^2{b}_{s^{\prime }t}^{\sigma }(n,\ell ).$$

$$B_{s-s}^{\sigma }(n,\ell )=N_{\ell s-s}^{\sigma }\frac{1}{N_{\ell s-s}^{\sigma }}{b}_{s-s}^{\sigma }(n\ell ),$$

$$=b_{s-s}^{\sigma }(n\ell ).$$

$$v(\theta ,\varphi ;T)=W(\theta ,\varphi )v_R(\theta ,\varphi ;T).$$

$$v_R({\theta }^{\prime },{\varphi }^{\prime };{\omega }^{\prime })=\sum _{s,t}{w}_{st}({\omega }^{\prime })Y_{st}({\theta }^{\prime },{\varphi }^{\prime }),$$

$${\theta }^{\prime }=\theta ,{\varphi }^{\prime }=\varphi +\Omega T,{\omega }^{\prime }=\omega -t\Omega .$$

$$v_R(\theta ,\varphi +\Omega T;T)=\frac{1}{2\pi }\int d{\omega }^{\prime }{e}^{i{\omega }^{\prime }T}\sum _s{w}_{st}({\omega }^{\prime })Y_{st}(\theta ,\varphi )e^{it\Omega T},$$

$$v_R(\theta ,\varphi ;T)=\frac{1}{2\pi }\int d\omega e^{i\omega T}\sum _s{w}_{st}(\omega -t\Omega )Y_{st}(\theta ,\varphi ).$$

$$v(\theta ,\varphi ;T)=W(\theta ,\varphi )\frac{1}{2\pi }\int d\omega e^{i\omega T}\sum _s{w}_{st}(\omega -t\Omega )Y_{st}(\theta ,\varphi ).$$

$$v_{st}(T)=L_{st}^{s^{\prime }{t}^{\prime }}\frac{1}{2\pi }\int d\omega e^{i\omega T}{w}_{s^{\prime }{t}^{\prime }}(\omega -t\Omega ),$$

Full text

helioseismology, waves — miscellaneous — catalogs — surveys

††journal: ApJ

1 Introduction

Rossby waves were first detected as large-scale weather patterns in Earth’s atmosphere (Rossby, 1939) and subsequently in the ocean (Chelton & Schlax, 1996). Theoretical analysis (Papaloizou & Pringle, 1978; Unno et al., 1989; Provost et al., 1981; Smeyers et al., 1981; Saio, 1982) suggested that rotating spherical fluids should show the presence of Rossby waves (commonly known as r-modes oscillations in the astrophysical context), in which the Coriolis force acts as a restoring force. The Sun, a differentially rotating spherical fluid, satisfies this condition. Long-term high-quality observations have encouraged several attempts in the past to detect Rossby waves (Kuhn et al., 2000; Ulrich, 2001; Sturrock et al., 2015). In all these earlier studies, the dispersion relation was not measured, a critical quantity in the identification of the governing physics. Recently, Löptien et al. (2018), using granulation tracking and the helioseismic technique of ring diagrams (Hill, 1988) measured a Rossby-dispersion relation, thereby unambiguously detecting sectoral modes (where the azimuthal order of the mode is the same as its harmonic degree) in the subsurface through analyses of six years of data from the Helioseismic and Magnetic Imager (HMI: Schou et al., 2012) onboard the Solar Dynamics Observatory. Subsequently, Liang et al. (2018) and Hanasoge & Mandal (2019), using two different methods, time-distance helioseismology (Duvall et al., 1993b) and normal mode coupling respectively, further validated the original detection. All the studies mentioned above have focused on Rossby waves in the solar surface and near-surface regions. Since Rossby waves exist in spherical rotating fluids, it can also in principle exist in solar atmosphere. McIntosh et al. (2017) using the tracking of coronal bright points has found evidence of Westward propagating large-scale Rossby-like waves in the solar corona. It may be useful to characterize the connection between solar internal Rossby waves and the measurements of McIntosh et al. (2017). One important difference between them is that solar internal Rossby waves are retrograde waves whereas Rossby like waves detected in solar corona by McIntosh et al. (2017) are prograde propagating. Zaqarashvili et al. (2010a, b, 2015); Dikpati et al. (2018a, b) studied magnetized Rossby waves in a magneto-hydrodynamic shallow water model of the solar tachocline in order to understand quasi-biennial oscillations and Rieger-type periodicities 111Solar activity follows a dominant $11$-year cycle. Two other timescales have also been observed in solar activity indicators, the first presenting timescales of around $150$-$160$ days - commonly known as Rieger-type periodicity and the second following an approximately-$2$ year periodicity, known as quasi-biennial oscillations. in solar activity. It is therefore important to characterize the properties of internal solar Rossby waves from observations in order to determine how it affects solar internal dynamics.

Normal-mode coupling, a seismic technique with an illustrious history in geophysics (see, e.g. Dahlen & Tromp, 1998), has seen limited adoption in helioseismology. Woodard (1989) first described distortion of eigenfunctions of the Sun due to latitudinal differential rotation. This method has subsequently been used by several authors e.g. Lavely & Ritzwoller (1992); Woodard (2006); Roth & Stix (2008); Schad et al. (2011, 2013); Woodard et al. (2013); Woodard (2014, 2016); Hanasoge (2017, 2018) in various studies ranging from meridional circulation to convection. Roth & Stix (2003); Woodard (2016); Hanasoge (2018); Hanasoge et al. (2017) have showed how this method may be used to glean information about time-varying, non-axisymmetric features in the Sun. The fundamental measurement comprises cross-correlated Fourier coefficients of the wavefield, $\langle\phi^{*}_{\ell^{\prime}m^{\prime}}(\omega^{\prime})\phi_{\ell m}(\omega)\rangle$, where angular brackets denote ensemble averaging, $\phi$ the line-of-sight Doppler velocity, $\ell$ and $m$ harmonic degree and azimuthal order respectively, and $\omega$ the temporal frequency. In the absence of perturbations, the expected value of this correlation is precisely obtained by considering leakage and by modelling modes as independently and identically distributed 222If an outcome of a random sampling is independent of the random variables that came before it, then we call the process independent and identically distributed in short i.i.d random processes. Depending on the structure of the perturbations, specific spatio-temporal wavenumbers show enhanced power above the background systematic floor. This constitutes a “detection” of a perturbation at those spatial and temporal scales. The difference between temporal frequencies $\sigma=\omega^{\prime}-\omega$ captures information about the perturbation time scale, whereas the difference between harmonic degrees, $\ell^{\prime}-\ell$, and azimuthal number, $m^{\prime}-m$ carry information about the length scale, and toroidal and poloidal nature of flows. Hanasoge & Mandal (2019) used this formalism to detect Rossby waves in the Sun with the first $2$ years of HMI data, thereby validating the methodology. In this work, we have extended the analysis of Hanasoge & Mandal (2019) with $12$ and $8$ years of solar oscillation data from Michelson Doppler Imager (MDI) onboard the Solar and Heliospheric Observatory and HMI respectively to estimate frequencies, line-widths and amplitudes of Rossby modes and compare them with parameters obtained from earlier studies by Löptien et al. (2018); Liang et al. (2018). Because we do not observe the entire Sun, spatial and temporal leakage affects our measurements. We perform several synthetic tests in order to characterize the effect of leakage, in order to be able to place faith in our inferences, e.g. to recover the depth profile of Rossby modes accurately from a synthetic measurement with added noise. This will help us further determine depth profiles of Rossby waves in the convection zone, which are still unknown and are active areas of research.

2 Data Analysis

The mode-coupling measurement procedure for Rossby waves is described in detail in Hanasoge & Mandal (2019) and Hanasoge (2018). We follow same notations as Hanasoge (2018) for convenience. The raw data are global time series of line-of-sight Doppler velocity projected on to spherical harmonics, i.e. $\phi_{\ell m}$ from MDI and HMI. We calculate the cross-correlation function $\langle\phi^{*}_{\ell m}(\omega)\phi_{\ell m+t}(\omega+\sigma)\rangle$. Analyzing this quantity at each frequency, $\omega$, is less tractable, and we therefore define B-coefficients as a linear least-square fit to the raw wavefield correlations,

[TABLE]

where $H$ is a weight function defined in Hanasoge & Mandal (2019) and Hanasoge (2018)

[TABLE]

where the first term on the right hand side of Equation (2) is a Wigner-$3j$ symbol and $L^{\ell^{\prime}m^{\prime}}_{\ell m}$ is the leakage matrix, which describes how spatial windowing in the data, i.e. arising from our limited vantage of the Sun, causes “leakage” from mode $(\ell,m)$ to another mode $(\ell^{\prime},m^{\prime})$. $B$-coefficients are calculated for all identified radial orders, $(n)$ and harmonic degrees in the range $\ell\in[10,180]$. $R^{\omega}_{\ell m}$ describes the power spectrum of a mode (labelled using three quantum numbers $n,\ell,m$) (Anderson et al., 1990; Duvall et al., 1993a),

[TABLE]

where $\omega_{n\ell m}$ is an eigenfrequency and $\Gamma_{n\ell}$ is the full width at half maximum. We use observed values, i.e. obtained through fits of spectra by the global-mode pipeline, for these parameters. We assume that errors in measurements of $\langle\phi^{*}_{\ell m}(\omega)\phi_{\ell m+t}(\omega+\sigma)\rangle$ are uncorrelated, have equal variances and zero expectation value. This allows us to invoke Gauss-Markov theorem which ensures that B-coefficients estimated from Equation 1 as an unregularized least square estimate of the cross correlation measurement are unbiased. Differential rotation advects features on the Sun at a variety of different speeds. In order to appropriately follow these perturbations, we adopt a co-rotating frame, i.e. apply a tracking rate. The temporal frequencies of the perturbation in the co-rotating frame, $\sigma$, will be transformed to $\sigma+t\Omega$ in the inertial frame, where $\Omega$ is rotation frequency of the Sun. A variety of tracking rates may be chosen since the Sun is differentially rotating: here, we choose $453$ nHz, which is the value at the equator. We extend the frequency range from our earlier work (Hanasoge & Mandal, 2019) to $\sigma\in[0,200]$ nHz in order to appreciate the spectra of the perturbations better.

3 Inversion

We assume Rossby waves are sufficiently well described by a toroidal flow, allowing it to be expressed as following

[TABLE]

where $(r,\theta,\phi)$ and $(\hat{\mathbf{r}},\hat{\bm{\theta}},\hat{\bm{\phi}})$ are radius, co-latitude, longitude and corresponding unit vectors respectively, $\bm{\nabla}_{h}$ is the horizontal covariant derivative and $Y_{st}$ is a spherical harmonic of degree $s$ and azimuthal-order $t$. $w^{\sigma}_{st}(r)$ determines the depth variation of Rossby waves.

We observe only the near side of the Sun, the angular extent that appears in the field of view of the telescope. This results in spatial and temporal leakage in the measurements (see section 4.1). B-coefficients as estimated in section (2) will encounter leakage from neighbouring harmonic degrees and azimuthal number (see Equation 26 of Hanasoge, 2018)

[TABLE]

where

[TABLE]

and

[TABLE]

$B$ and $b$-coefficients refer to different quantities and are only in special cases (as discussed below) identical. $B^{\sigma}_{st}$ and $b^{\sigma}_{st}$ are also functions of radial order, $n$. In this work we consider the coupling between modes with same radial order and therefore omit it in the expression for notational brevity. $\mathcal{K}_{n\ell}$ is the sensitivity kernel for the mode $(n,\ell)$ and $f_{\ell^{\prime}-\ell,s}$ is obtained from an asymptotic analysis of the kernels (Vorontsov, 2011; Hanasoge, 2018; Woodard, 2014). This asymptotic kernel is only valid when $s\ll\ell$ or $s\ll\ell^{\prime}$ (see for example Brussaard & Tolhoek, 1957; Vorontsov, 2011). Hanasoge (2018) have shown a comparison between kernels obtained from asymptotic analysis and the exact expression (Figure 10 of that paper) and found that the comparison gets slightly worse when $s\approx\ell$ but is otherwise very accurate when $s<\ell$. In this work, we search for coupling between p-modes with harmonic degree in the range $\ell\in[10,180]$, and are interested in Rossby modes with harmonic degree $s\leq 20$, thereby justifying the use of asymptotic kernels. In case of full-sphere observations, we show in appendix A that the measured B-coefficients from Equation (1) will be reduced to Equation (7), a simple relation with which to invert for the velocity profile $w^{\sigma}_{st}$.

Inverting Equation (5) is complicated in the general case, as explained in Hanasoge (2018) due to leakage contributions from neighbouring modes into the target measurement channel. For solar Rossby waves, we know from earlier studies by Löptien et al. (2018); Liang et al. (2018); Hanasoge & Mandal (2019) that non-sectoral modes $(s\neq|t|)$ are either absent or appear at amplitude levels well below the detection threshold. Due to this condition, Equation 5 would be simplified after substituting $b^{\sigma}_{st}\approx\delta_{s,-t}b^{\sigma}_{s-s}$

[TABLE]

Although Equation 8 is much simpler than Equation 5, we still need to show that we can separate out desired specific spatial scales using our measurement. We follow two approaches as discussed below for that. In our first approach, we do not consider leakage and invert Equation (7), assuming $b^{\sigma}_{s,-s}$ is the same as $B^{\sigma}_{s,-s}$, as considered in our earlier work (Hanasoge & Mandal, 2019). This approximation will only hold if there is no leakage, i.e., $L^{\ell^{\prime}m^{\prime}}_{\ell m}=\delta_{\ell\ell^{\prime}}\delta_{mm^{\prime}}$. The approximation affects our inferences and to quantify this, we perform synthetic tests. In our second approach we consider leakage and a discussion of which is presented in section  (4.2). We choose two inversion techniques, Optimally Localized Averaging, (OLA; Backus & Gilbert, 1967) and Regularized Least Square, (RLS) for this work.

3.1 OLA

In OLA, the inverted flow profile at depth $r_{0}$ is written as a linear combination of all the measurements as following,

[TABLE]

where $\alpha$ needs to be determined by minimizing the following misfit function

[TABLE]

$\lambda$ is the regularization parameter. $\sum_{n,\ell}\alpha_{n\ell;r_{0}}f_{0,s}K_{n\ell}$ determines the averaging kernel at depth $r_{0}$. Weight function $(r-r_{0})^{2}$ in the first term of Equation (10) ensures that the averaging kernel obtained after minimizing Equation (10) is large valued at depth $r_{0}$ and small elsewhere. We can choose to perform the inversion separately for each frequency bin or once for all frequency bins. In this technique, we choose the latter approach: ${\mathcal{N}}_{n\ell}$ is obtained from $\mathcal{N}^{\sigma}_{n\ell}$ after summing over all the frequency bins, $\sigma$, where $\mathcal{N}^{\sigma}_{n\ell}$ is the diagonal component of the noise covariance matrix. In the next section we perform an inversion for each frequency bin separately using the RLS method.

3.2 RLS

In this method, we expand the flow profile in the B-spline basis

[TABLE]

where $B_{k}$ is the B-spline basis function of order three. We choose total $50$ knots for the inversion up to depth $0.1R_{\odot}$. We determine $\beta$ by minimizing the following misfit function

[TABLE]

where we consider second-derivative smoothing where $\lambda$ is the regularization parameter. The Equation (12) may be minimized if we solve the following system of equations

[TABLE]

[TABLE]

where $r_{k}$ denotes points on the radial grid. As opposed to earlier methods, we must perform inversions for separate frequency bins in RLS.

4 Results

4.1 Modelling leakage

Hanasoge & Mandal (2019) considered a frequency analysis $\sigma$ in the range $(0.0,0.5)\,\mu$Hz. In this work, we expand the range of $\sigma$ up to $2\,\mu$Hz and show in Figure (1) that Rossby modes leak into higher frequencies, much as observed by Liang et al. (2018). The sectoral mode of harmonic degree $s$ and corresponding frequency $\sigma_{s}=2\Omega/(s+1)$ leaks into a harmonic degree $s+2$ at frequency $\sigma_{s}+2\Omega$. Spatial windowing of the rotating Sun results in simultaneous spatial and temporal leakage, which we verify through a synthetic test.

We consider sectoral modes of Rossby waves for odd harmonic degrees. Rossby waves in the rotating Sun are multiplied by a spatial window function - which is unity for the visible portion of the disk, and otherwise zero. We perform spherical harmonic and temporal Fourier transforms in order to estimate $v_{s}(\omega)$. Since we only detect sectoral retrograde modes, i.e. $t=-s$, we use both $t$ and $s$ equivalently to denote Rossby waves. Analytical calculations suggest (see Equation B12 in Appendix B) that the observed $v_{s}(\omega)$ contains contributions from neighbouring modes $s^{\prime}$ due to leakage $L_{s}^{s^{\prime}}$,

[TABLE]

The Delta function is invoked assuming that Rossby modes have power close to the classical dispersion frequency

[TABLE]

where $\sigma_{t}$ is the frequency of the Rossby mode with azimuthal number $t$ in a co-rotating frame with rotation frequency $\Omega$. From Equation (15) we see that neighbouring modes $s^{\prime}$ will corrupt the measurements of our desired mode $s$ by an amount that depends on the value of leakage $L_{s}^{s^{\prime}}$.

As explained in Appendix B, the dispersion relation of Rossby waves either leaks into $\Omega$ or $2\Omega$ depending on the value of $L_{s}^{s^{\prime}}$. In order to appreciate which of $L_{s}^{s+1},$ or $L_{s}^{s+2}$ is more significant, we calculate $B$-coefficients using Equation (8) without considering tracking and obtain the quantity $(\sum_{n}(-1)^{\ell}B^{\sigma}_{st}(n,\ell))^{2}$. From the left panel in Figure (2), we see that power from odd $s$ leaks to odd $s$ and not to even $s$ for our measurements, which implies that $L_{s}^{s+2}$ is more significant than $L_{s}^{s+1}$. Because of this reason, we observe leakage at $\omega_{t}+2\Omega$ but not at $\omega_{t}+\Omega$. The right panel of Figure (2) shows the effect of leakage if we track the data, similar to the spectrum observed in Figure (1). There might be other systematics, e.g., P or B-angle corrections, that can affect the properties of leakage as discussed above.

4.2 Synthetic inversions taking into account leakage

Hanasoge (2018) have shown that leakage complicates the inversion using Equation (5). Therefore in our earlier work (Hanasoge & Mandal, 2019), we performed inversions using Equation (7) for simplicity’s sake instead of the full Equation (5). This simplifying assumption might diminish the accuracy in retrieving the depth profiles of Rossby waves. In our synthetic test, we use cubic polynomials to characterize the depth profiles, with the condition that they are set to zero at depths $0.9R_{\odot}$ and below. We subsequently calculate the $B$-coefficient using Equation (5). In our first approach, we ignore leakage and assume $b^{\sigma}_{st}$ is same as $B^{\sigma}_{st}$ and then use Equation (7) to invert for the profile $w^{\sigma}_{st}$.

Next, we take leakage into account and proceed with the following approach. We write Equation (8) using the following compressed form

[TABLE]

where

[TABLE]

is the new sensitivity kernel which relates observed B-coefficients to properties of Rossby waves $w^{\sigma}_{s\,-s}$. The off-diagonal terms, $\Theta_{s}^{s^{\prime}}$ quantify contributions from neighbouring modes $(s^{\prime},-s^{\prime})$ to our desired mode $(s,-s)$ in the observed B-coefficient. The term $\Theta_{s}^{s^{\prime}}(\sigma)$ depends on the leakage matrix $L_{\ell m}^{\ell^{\prime}m^{\prime}}$. We have shown in section (4.1) that leakage from neighbouring modes $s^{\prime}$ occurs at the same temporal frequency as that of the mode $\sigma_{s^{\prime}}$ if we do not employ tracking. If the data are tracked, this leakage moves to higher temporal frequencies. Since we are performing inversions at each frequency bin in the range $[0.0,0.5]\,\mu$Hz and since the contribution to our desired frequency bins from neighbouring modes is negligible in that range, we consider diagonal terms only, i.e. $\Theta_{s}^{s}$ in Equation (18) for the inversion.

[TABLE]

We apply Equation (19) to invert for $B$-coefficients, estimated using Equation  (5). In the left panel of Figure (3), we compare our inversion results with the original profile, demonstrating that choosing diagonal terms in Equation (19) improves the inverted profile. This test also shows that choosing Equation (7) is not a particularly bad assumption for this problem, and that our inferred amplitude might at worst be off by a factor of a few. It is the kernel that is responsible for differences in the inferred amplitudes. To compare the two kernels, we plot $f_{0,s}\mathcal{K}_{n\ell}$ and $\Theta_{s}^{s}$ in the right panel of Figure (3) for a mode with radial order $n=2$ and harmonic degree, $\ell=131$. It can be seen in Figure (3) that though the shape of the kernels remains same in two cases, the values are different, resulting in changes in the inferred amplitude. In order to test the inversion algorithm with noise added to the measurement, we choose an artificial profile that goes to zero at the base of the convection zone. The magnitude at the surface is set to $4$ ms*-1*, which is close to the observed value. We calculate $B$-coefficients in a similar manner as described above using Equation (5) and add random Gaussian noise in proportion to the observed level, ${\mathcal{N}}_{n\ell}^{\sigma}$. We then perform inversions assuming Equation (7) and compare our inversion results with the original profile in Figure (4). In order to get error in the inverted profile we repeat this process multiple times and estimate the standard deviation of these values, thereby obtaining the error in the inferential profile.

4.3 Frequencies and line-widths of r-modes

After validating our method we estimate B-coefficients as described in section (2) from the observed oscillation data. We divide $12$ years of MDI data into three four-years chunks and 8 years of HMI data into two four-years chunks. We analyze each chunk of data sets separately and do the inversion to obtain, $w^{\sigma}_{st}$ from the measured B-coefficients. We then take average of these results and plot the average power spectra which is shown in Figure (5). In order to quantify frequencies, line-widths and amplitude of these modes, we fit a Lorentzian function plus a constant background

[TABLE]

to $|w_{st}|^{2}$ for each $s$. Here $A$ is the maximum amplitude of the Lorentzian, $\omega_{0}$ is the central frequency, $\Gamma$ is the full width at half maximum, $B$ is constant background power. We use the curve_fit module implemented in scipy.optimize to fit the power spectrum. We have tabulated values of these parameters for all modes starting from $s=1$ to $s=15$ obtained through analyses of HMI and MDI data in Tables (1) and (2) respectively. Fitted spectrum for HMI and MDI are shown in Figure (6) and (7) respectively. The $s=3$ mode parameters obtained from HMI data are similar to the values reported by Löptien et al. (2018), whereas the MDI analysis is similar to the findings of Liang et al. (2018). The $s=13$ mode obtained from the analysis of HMI data does not clearly stand out and so we do not fit this mode. Amplitudes of all the modes are varying with harmonic degrees and becoming very small for harmonic degrees $s\geq 13$. Since our analysis period covers cycles $23$ and $24$, differences in parameter values between Tables (1) and (2) (obtained from analysis of HMI and MDI data respectively) may in principle carry information about the solar cycle dependence of Rossby modes. However, this interpretation is complicated because of our use of data from two different instruments that we have not cross calibrated. Liang et al. (2017) have performed this in the context of time-distance helioseismology. A similar task is required for normal mode coupling in order to combine MDI and HMI analyses.

4.4 Is the mode of degree $s=1$ due to systematics in the method?

Hanasoge & Mandal (2019) detect the $s=1$ Rossby mode. Since the tracking rate is same as the frequency of the mode, $s=1$, systematics in our method might induce spurious power in the related spatial and temporal frequency bins. To investigate if it is due to tracking, we choose different values for it and calculate the B-coefficient for each case. If the tracking rate changes, so will the frequency of each mode in the power spectrum plot in Figure ( 5). We plot the frequency of the mode, $s=1$ versus tracking rate and compare it with the theoretical value in Figure (8). We also see that the $s=1$ mode is leaking to mode $s=3$ in Figure (1), which is unlikely to occur if it were to be due to systematics in the method. These arguments are in favour of the theory that the $s=1$ power is associated with Rossby waves and not due to systematics. We also plot the power for harmonic degree, $s=1$ and azimuthal number, $t=1$ in the right panel of Figure (8). We do not find extra power close to the frequency $453$ nHz. Nevertheless, caution should be taken when interpreting this mode as we might not have accounted for all the systematics in our method. One possibility might be centre-to-limb systematic which is viewed as an outflow from the disk centre to the limb. This is a poloidal flow. Nevertheless, if there is a component of this poloidal outflow that leaks to a toroidal component with harmonic degree, $s=1$ and azimuthal number, $t=-1$, spurious power at this temporal and spatial scale would show up in the spectrum. This fictitious flow is static from an inertial frame and therefore in a co-rotating frame its frequency will be equal to the tracking rate $\Omega$, which is same as the Rossby mode, $(s,t)=(1,-1)$.

5 Conclusion

We have extended the work by Hanasoge & Mandal (2019) by analyzing $8$ years of HMI data by considering all observed oscillation modes of harmonic degree in the range $[10,180]$ to measure frequencies, line-widths and amplitude of Rossby modes. We also analyze $12$ years of SOHO/MDI data from years, $1999$ to $2011$. In this work, we do not combine these two data sets, since a more sophisticated analysis to correct for P and B-angle variations is required (as demonstrated by Liang et al., 2017, in the context of time-distance helioseismology). Therefore, we show results from the two data sets separately.

In order to understand how Rossby mode eigenfunctions vary with depth, we first validate our inversion technique by recovering a synthetic profile that we introduce. We first use Equation (7) and perform the inversion. Next we consider the diagonal term in Equation (17) for the inversion. In Figure (3), we compare the inverted and original profiles. Choosing the diagonal term for the inversion provides a higher quality inference than when simply choosing Equation (7).

We performed tests to verify if the $s=1$ mode that we detect is due to systematics in our measurement since the frequency of this mode is same as tracking rate we choose, $\Omega$. We choose different tracking rates and find that the observed frequency of the $s=1$ mode follows the theoretical prediction for that tracking rate. We also see that $s=1$ leaks to $s=3$ (see Figure 1), which is unlikely to occur if it were to be due to systematics in our method. However we remain skeptical and will not declare detection of the $s=1$ mode.

Our approach of modelling leakage and measurement technique differs from that of Schad et al. (2013), who studied meridional circulation in the Sun using a variant of the measurement we apply here. Here, we consider coupling of identical-harmonic-degree p-modes, which limits us to the study of Rossby waves of odd harmonic degrees. In order to investigate the even harmonic degrees, we need to consider coupling between p-modes of different harmonic degrees, with $\delta\ell=1,3,\ldots$ (as $\delta\ell$ increases, frequency separation also increases, which in turn decreases the sensitivity of the coupling, as may be seen from Equations 1 and 2). Additionally, in order to image meridional circulation and differential rotation, we have to estimate $B$-coefficients as defined in Equation 1 at $\sigma=0$ (since we assume that these features only weakly evolve in time). We also mention in section 4.2 how the particular problem of Rossby waves in which we are interested becomes simplified because of its finite $\sigma$. Our future attempts will be to address more general problems in helioseismology using this technique.

Appendix A Partial to full sphere observation by ignoring leakage

In case of full-sphere observations of the Sun, there is no leakage, i.e. $L_{\ell m}^{\ell^{\prime}m^{\prime}}=\delta_{\ell\ell^{\prime}}\delta_{mm^{\prime}}$. The general expression for the $B$-coefficient is

[TABLE]

If we substitute $L_{\ell m}^{\ell^{\prime}m^{\prime}}=\delta_{\ell\ell^{\prime}}\delta_{mm^{\prime}}$ into Equation A1 we obtain

[TABLE]

Note that we have put dependence on radial order $n$ back in the last step of the above equation. We currently observe only sectoral Rossby modes, i.e. $s=-t$ in the Sun. If we apply this condition to the above equation and using Equation 6, we arrive at

[TABLE]

Note that, to arrive at the final Equation A3, we have assumed that only sectoral modes of Rossby waves exist.

Appendix B Spatial and temporal leakage

We observe only a part of the solar disk that comes into the field of view of the telescope. If $v_{R}(\theta,\phi;T)$ is the velocity field due to Rossby waves, our measurement will be $v(\theta,\phi;T)$ and these two fields are related by the following expression

[TABLE]

$W(\theta,\phi)$ is the window function unity over the visible part of the solar disk and otherwise zero. In our convention, $T$ denotes time and $t$ is used for azimuthal number. We assume that the velocity field for Rossby waves, $v_{R}(\theta,\phi;T)$, is scalar for simplicity. The velocity field of Rossby waves in the co-rotating frame may be written as

[TABLE]

where $Y_{st}$ is the spherical harmonic with azimuthal number $t$ and harmonic degree $s$. $w_{st}(\omega)$ is the inferred toroidal velocity field at that spatio-temporal frequency bin and radius. The prime denotes the coordinates of a co-rotating frame at angular frequency $\Omega$ rotating (with respect to the inertial frame). Coordinate and frequency transformations from one to other coordinate system are given by

[TABLE]

We omit the prime on $\theta$ in all equations below. Substituting Equation (B3) into (B2) and using (B1), we obtain

[TABLE]

Using Equations (B4) and (B1), we have

[TABLE]

Performing a spherical-harmonic transform of Equation (B5), we obtain

[TABLE]

where $L_{st}^{s^{\prime}t^{\prime}}$ denotes leakage from mode $(s,t)$ to $(s^{\prime},t^{\prime})$,

[TABLE]

Access to the full-sphere observation would have meant $L_{st}^{s^{\prime}t^{\prime}}=\delta_{ss^{\prime}}\delta_{tt^{\prime}}$. After performing a temporal Fourier transform of Equation (B6), we obtain

[TABLE]

Since we only detect sectoral modes of Rossby waves, we drop $s$ or $t$ keeping in mind that $t=-s$. With this substitution, Equation (B8) becomes

[TABLE]

The frequencies of Rossby waves in the inertial frame follow the relation

[TABLE]

where $\sigma_{t}$ is the frequency of Rossby wave modes with azimuthal order $t$ which obeys Equation (16). The power of these modes lies close to these frequencies (since the linewidths are small, see, e.g., Table 1). For simplicity, we choose a delta function

[TABLE]

After substituting Equation (B11) into (B9), we obtain

[TABLE]

At frequencies close to $\sigma_{t}+t\Omega$ we have

[TABLE]

which is non zero only when $t=t^{\prime}$

[TABLE]

which in turn is the classical Rossby-wave dispersion relation in an inertial frame. For frequencies $\sigma_{t+1}+(t+1)\Omega$ and $\sigma_{t+2}+(t+2)\Omega$, we recover

[TABLE]

We would see power at frequency $\sigma_{t+1}+\Omega$ in a co-rotating frame with azimuthal order $t$ if $L_{s}^{s+1}$ were to be significant. Similarly, power at $\sigma_{t+2}+2\Omega$ is observed if $L_{s}^{s+2}$ were to be large.