TY - JOUR
AU - Erb,, Wolfgang
AB - Abstract For sampling values along spherical Lissajous curves we establish a spectral interpolation and quadrature scheme on the sphere. We provide a mathematical analysis of spherical Lissajous curves and study the characteristic properties of their intersection points. Based on a discrete orthogonality structure we are able to prove the unisolvence of the interpolation problem. As basis functions for the interpolation space we use a parity-modified double Fourier basis on the sphere that allows us to implement the interpolation scheme in an efficient way. We further show that the numerical condition number of the interpolation scheme displays a logarithmic growth. As an application, we use the developed interpolation algorithm to estimate the rotation of an object based on measurements at the spherical Lissajous nodes. 1. Introduction In magnetic resonance imaging, motions of the scanned subject during the imaging process cause artifacts in the reconstruction. One concept to detect and correct subject motions in three dimensions is based on the additional measurements of spherical navigator echoes (Welch et al., 2002; Costa et al., 2010). These measurements are performed along sampling trajectories on a spherical shell and used to estimate rotations and translations of the scanned subject. Particular promising trajectories for such navigator measurements are spherical Lissajous curves (Ullisch, 2012). A spherical Lissajous curve is given in parametric form as \begin{equation} \boldsymbol{\ell}^{(\boldsymbol{m})}_{\alpha}(t) = \Big( \sin(m_2 t) \cos (m_1 t - \alpha \pi), \ \sin(m_2 t) \sin (m_1 t - \textstyle \alpha \pi), \ \cos(m_2 t) \Big), \quad t \in \mathbb{R}, \end{equation} (1.1) with a frequency vector |$\boldsymbol{m} = (m_1, m_2) \in{{\mathbb{N}}}^{2}$| and a rotation parameter |$\alpha \in{{\mathbb{R}}}$|⁠. The curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| lies in the unit sphere |${\mathbb{S}}^2 = \{\boldsymbol{x} \in{{\mathbb{R}}}^3:\; x_1^2+x_2^2+x_3^2 = 1\}$| of the three-dimensional space |${{\mathbb{R}}}^3$|⁠. Similar as for bivariate Lissajous curves (Erb, 2016; Erb et al., 2016; Dencker & Erb, 2017a,b), the curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| describes a superposition of a latitudinal and a longitudinal harmonic motion determined by the frequencies |$m_1$| and |$m_2$|⁠. The goal of this article is to provide a mathematical analysis of spherical Lissajous curves and to study their role as generating curves for spherical interpolation. Of particular interest for our analysis are the intersection points |$\textbf{LS}^{(\boldsymbol{m})}$| of one or more spherical Lissajous curves. These Lissajous nodes provide a good measure of how densely the curves cover the sphere. Further, these nodes are relevant for applications. In Ullisch (2012), the intersection points of spherical Lissajous curves are used to correct spin–spin relaxation effects for the navigator measurements. In the upcoming sections, we will give two different characterizations of the Lissajous nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. Particularly interesting for us and for applications is the case when the two frequencies |$m_1$| and |$m_2$| are relatively prime and |$m_2$| is even. In this case, the nodes |$\textbf{LS}^{(\boldsymbol{m})}$| can be described as time equidistant samples along a single spherical Lissajous curve. The restriction to even numbers |$m_2$| guarantees that the two poles of |${\mathbb{S}}^2$| are included in |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. The second characterization in Section 3 is based on particularly defined index sets |$\textbf{I}^{(\boldsymbol{m})}$| and allows an explicit description of the nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. It includes also the case when |$m_1$| and |$m_2$| are not relatively prime and when more than one Lissajous curve is needed to generate the nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. The main task of this manuscript is to use the spherical Lissajous nodes |$\textbf{LS}^{(\boldsymbol{m})}$| to derive a novel scheme for interpolation and quadrature on the sphere. To this end, we transfer concepts and techniques developed in Dencker & Erb (2017a,b) for multivariate polynomial interpolation on Lissajous–Chebyshev nodes in the hypercube |$[-1,1]^{\mathsf{d}}$| to a corresponding setting in spherical coordinates. As for the multivariate Lissajous–Chebyshev points, a main step in the proof of the spherical interpolation scheme is a discrete orthogonality structure linked to the spherical Lissajous nodes. This discrete orthogonality structure will be derived in Section 4. As a basis system for the interpolation on the Lissajous nodes we will use a parity-modified double Fourier basis in spherical coordinates. These basis functions were introduced in the 70’s (Merilees, 1973; Orszag, 1974) as a stable alternative to the spherical harmonics and the Robert functions. Since then, they were used in a series of applications on the sphere as for instance described in Boyd (2000, Section 18.27) and Boer & Steinberg (1975), Boyd (1978), Fornberg (1995), Ganesh et al. (1998) and Townsend et al. (2017). As these functions are directly built on a Fourier series they are very well suited for computational purposes. Compared to spherical harmonics there are however some issues at the poles of the sphere that have to be treated properly. For a more detailed discussion on different aspects of this basis system we refer to the treatises given in Boyd (1978, 2000). For the actual work the concrete form of the parity-modified Fourier basis plays a crucial role. We will establish a close link between interpolation on the nodes of the spherical Lissajous curves and this basis system. This discussion will lead to Theorem 5.2 in which we prove the uniqueness of the interpolation in spaces spanned by the double Fourier basis. In Section 6, we will discover that the mentioned structure of the basis functions leads to an efficient implementation of the interpolation scheme in terms of a double Fourier transform. In Section 7 we will further give a short description of the numerical condition number and the convergence of the interpolation scheme. We will see that, similar to spectral methods on the hypercube |$[-1,1]^{\mathsf{d}}$|⁠, the numerical condition displays a slow logarithmic growth and that the interpolant converges fast if the data values are derived from smooth functions. Finally, we will provide some numerical experiments and present an idea on how this novel interpolation scheme on the sphere can be applied to estimate rotations of a three-dimensional object based on measurements along spherical Lissajous curves. 2. Spherical Lissajous curves In the first step, we want to derive some fundamental properties of the spherical Lissajous curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| defined in (1.1). In particular, we are interested in its minimum period and in the number and type of its self-intersection points. Two examples of spherical Lissajous curves with their intersection points are illustrated in Fig. 1. Fig. 1. View largeDownload slide Illustration of two Lissajous curves and its intersection points described in Proposition 2.1. The intersection points of the curve are marked with black dots. The circled dots describe the two poles of |${\mathbb{S}}^2$|⁠. Fig. 1. View largeDownload slide Illustration of two Lissajous curves and its intersection points described in Proposition 2.1. The intersection points of the curve are marked with black dots. The circled dots describe the two poles of |${\mathbb{S}}^2$|⁠. If the frequencies |$m_1$| and |$m_2$| of |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| are relatively prime, Proposition 2.1 below implies that the minimum period of |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| is |$2\pi $|⁠. In general, if |$g = \textrm{gcd}(\boldsymbol{m})$| denotes the greatest common divisor of |$m_1$| and |$m_2$|⁠, then |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| can be rewritten as |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t) = \boldsymbol{\ell }^{(\boldsymbol{m}/g)}_{\alpha }( g t)$| and the minimum period of |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| is |$2 \pi /g$|⁠. For the description of spherical Lissajous curves it is therefore enough to restrict ourselves for the moment to tuples |$\boldsymbol{m}$| of relatively prime numbers. To extract the self-intersection points of the curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| we consider for |$t \in [0,2\pi )$| the sets |$\mathcal{A}^{(\boldsymbol{m})}(t) = \{ s \in [0,2\pi ): \ \boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(s) = \boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t) \}$| and the sampling points \begin{align} t^{(\boldsymbol{m})}_{l} &= \frac{l \pi}{m_1 m_2}, \quad l \in \{0,1, \ldots, 2m_1m_2-1\},\qquad\!\!\! \end{align} (2.1) \begin{align} t^{(\boldsymbol{m})}_{l+\frac12} &= \frac{\big(l + \frac{1}{2}\big) \pi}{m_1 m_2}, \quad l \in \{0,1, \ldots, 2m_1m_2-1\}. \ \end{align} (2.2) Proposition 2.1 Let |$\textrm{gcd}(\boldsymbol{m}) = 1$|⁠. If |$m_2$| is even, then \begin{equation*} \begin{array}{lll} \textrm{(i)} & \# \mathcal{A}^{(\boldsymbol{m})}(t) = m_2 & \textrm{if}\quad t \in \{\,t^{(\boldsymbol{m})}_{l}\,| \ l\in \{0,\ldots,2m_1m_2-1\}, \ l \equiv 0 \mod m_1 \},\\ \textrm{(ii)} & \# \mathcal{A}^{(\boldsymbol{m})}(t) = 2 & \textrm{if}\quad t \in \{\,t^{(\boldsymbol{m})}_{l}\,| \ l\in \{0,\ldots,2m_1m_2-1\}, \ l \not\equiv 0 \mod m_1 \}, \\ \textrm{(iii)} & \# \mathcal{A}^{(\boldsymbol{m})}(t) = 1 & \textrm{if}\quad t \in [0,2\pi) \setminus \{\,t^{(\boldsymbol{m})}_{l}\,|\, l\in \{0,\ldots,2m_1m_2-1\}\}. \end{array} \end{equation*} If |$m_2$| is odd, then \begin{equation*} \begin{array}{lll} (\textrm{i})^{\prime} & \# \mathcal{A}^{(\boldsymbol{m})}(t) = m_2 & if \quad t \in \{\,t^{(\boldsymbol{m})}_{l}\,| \ l\in \{0,\ldots,2m_1m_2-1\}, \ l \equiv 0 \mod m_1 \}, \\ (\textrm{ii})^{\prime} & \# \mathcal{A}^{(\boldsymbol{m})}(t) = 2 & if \quad t \in \{\,t^{(\boldsymbol{m})}_{l+\frac{1}{2}}\,| \ l\in \{0,\ldots,2m_1m_2-1\} \}, \\ (\textrm{iii})^{\prime} & \# \mathcal{A}^{(\boldsymbol{m})}(t) = 1 & \textrm{for all other}\ t\ \in [0,2\pi). \end{array} \end{equation*} Remark 2.2 If |$m_1$| and |$m_2$| are relatively prime we obtain as a consequence of Proposition 2.1 that the minimum period of |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| is |$2 \pi $|⁠. The points |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t)$| with |$\# \mathcal{A}^{(\boldsymbol{m})}(t) = m_2$| correspond to the north or the south pole of the sphere. Hence, in every period |$[0,2\pi )$| the curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| traverses both poles |$m_2$| times. All other points |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t)$| with |$\# \mathcal{A}^{(\boldsymbol{m})}(t) = 2$| correspond to nonpolar double points of the curve on the sphere, i.e., they are traversed twice by the curve as |$t$| varies from |$0$| to |$2 \pi $|⁠. Depending on whether |$m_1$| is odd or even, we get a different number of self-intersection points for |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$|⁠. These numbers are summarized in Table 1. Table 1 Number and type of intersection points (IPs) for the curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| if |$m_1$|⁠, |$m_2$| are relatively prime. For general |$\boldsymbol{m}$|⁠, the corresponding numbers are obtained by considering the curve |$\boldsymbol{\ell }^{(\boldsymbol{m}/g)}_{\alpha }$| instead (⁠|$g = \textrm{gcd}(\boldsymbol{m})$|⁠) Curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| Number of IPs Type of IPs |$m_2$| even |$m_2 (m_1 -1) + 2$| $$\begin{array}{l} \text{two poles, traversed $m_2$ times in one period,} \\ \text{$m_2 (m_1 -1)$ nonpolar double points} \end{array}$$ |$m_2>1$| odd |$ m_1 m_2 + 2$| $$\begin{array}{l} \text{two poles, traversed $m_2$ times in one period,} \\ \text{$m_1 m_2$ nonpolar double points} \end{array}$$ |$m_2=1$| |$ m_1 $| $$\begin{array}{l} \text{$m_1$ nonpolar double points} \end{array}$$ Curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| Number of IPs Type of IPs |$m_2$| even |$m_2 (m_1 -1) + 2$| $$\begin{array}{l} \text{two poles, traversed $m_2$ times in one period,} \\ \text{$m_2 (m_1 -1)$ nonpolar double points} \end{array}$$ |$m_2>1$| odd |$ m_1 m_2 + 2$| $$\begin{array}{l} \text{two poles, traversed $m_2$ times in one period,} \\ \text{$m_1 m_2$ nonpolar double points} \end{array}$$ |$m_2=1$| |$ m_1 $| $$\begin{array}{l} \text{$m_1$ nonpolar double points} \end{array}$$ View Large Table 1 Number and type of intersection points (IPs) for the curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| if |$m_1$|⁠, |$m_2$| are relatively prime. For general |$\boldsymbol{m}$|⁠, the corresponding numbers are obtained by considering the curve |$\boldsymbol{\ell }^{(\boldsymbol{m}/g)}_{\alpha }$| instead (⁠|$g = \textrm{gcd}(\boldsymbol{m})$|⁠) Curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| Number of IPs Type of IPs |$m_2$| even |$m_2 (m_1 -1) + 2$| $$\begin{array}{l} \text{two poles, traversed $m_2$ times in one period,} \\ \text{$m_2 (m_1 -1)$ nonpolar double points} \end{array}$$ |$m_2>1$| odd |$ m_1 m_2 + 2$| $$\begin{array}{l} \text{two poles, traversed $m_2$ times in one period,} \\ \text{$m_1 m_2$ nonpolar double points} \end{array}$$ |$m_2=1$| |$ m_1 $| $$\begin{array}{l} \text{$m_1$ nonpolar double points} \end{array}$$ Curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| Number of IPs Type of IPs |$m_2$| even |$m_2 (m_1 -1) + 2$| $$\begin{array}{l} \text{two poles, traversed $m_2$ times in one period,} \\ \text{$m_2 (m_1 -1)$ nonpolar double points} \end{array}$$ |$m_2>1$| odd |$ m_1 m_2 + 2$| $$\begin{array}{l} \text{two poles, traversed $m_2$ times in one period,} \\ \text{$m_1 m_2$ nonpolar double points} \end{array}$$ |$m_2=1$| |$ m_1 $| $$\begin{array}{l} \text{$m_1$ nonpolar double points} \end{array}$$ View Large Proof. We use the equivalence relation |$s\eqsim t$| to denote that |$t-s\in 2\pi \mathbb{Z}$|⁠. We consider first all |$t \in [0,2\pi )$| such that |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t)$| is one of the poles |$(0,0,v)$|⁠, |$v\in \{-1,1\}$|⁠, of the unit sphere. By the definition (1.1) of the Lissajous curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$|⁠, the identity |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t)=(0,0, v)$| holds if and only if \begin{equation*} m_2 t \eqsim \frac{1-v}{2} \pi \quad \textrm{for}\ v \in \{-1,1\}, \end{equation*} i.e., if and only if |$t = t^{(\boldsymbol{m})}_{l}$| with |$l\in \{0,m_1, 2m_1, \ldots ,(2m_2-1)m_1\}$|⁠. Further, we have in this case |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t^{(\boldsymbol{m})}_{l})=(0,0, 1)$| if |$l/m_1$| is even and |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t^{(\boldsymbol{m})}_{l})=(0,0, -1)$| if |$l/m_1$| is odd. This yields the statements (i) and (i)′ of Proposition 2.1. We consider now the case that |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t)$| is not one of the poles. By the definition (1.1) of the curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t)$|⁠, we have |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(s) = \boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t)$| if and only if \begin{equation*}\cos(m_2 s) = \cos(m_2 t) \quad \textrm{and} \quad \sin(m_2 s) \left( \begin{array}{l} \cos(m_1 s - \alpha \pi) \\ \sin(m_1 s - \alpha \pi) \end{array} \right) = \sin(m_2 t) \left( \begin{array}{l} \cos(m_1 t - \alpha \pi) \\ \sin(m_1 t - \alpha \pi) \end{array} \right).\end{equation*} In the first formula we get equality precisely if |$m_2 s \eqsim v m_2 t$| for some |$v \in \{-1,1\}$|⁠. Plugging this relation into the second formula, we see that |$\sin (m_2 s) = v \sin (m_2 t)$| and, thus, that we get equality in the second formula exactly if |$m_1 s \eqsim m_1 t + \frac{1-v}{2} \pi $|⁠. In total, we can conclude that |$s\in \mathcal{A}(t)$| if and only if \begin{equation} m_2( s - v t) \eqsim 0 \quad \textrm{and} \quad m_1 ( s - t) + \frac{1-v}{2} \pi \eqsim 0 \qquad \textrm{for}\ v \in \{-1,1\} \end{equation} (2.3) is satisfied. Since |$m_1$| and |$m_2$| are relatively prime, Bézout’s lemma gives two integers |$a, b\in \mathbb{Z}$| such that |$a m_1 + b m_2 = 1$|⁠. The two conditions in (2.3) imply \begin{equation} s \eqsim t - b m_2 (1-v) t - a \frac{1-v}{2} \pi, \end{equation} (2.4) and, thus, |$s \eqsim t$| for |$v = 1$|⁠. For |$v = -1$|⁠, the two conditions in (2.3) imply \begin{align*} 2 m_1 m_2 s\eqsim m_1 m_2 (s + t) + m_2 m_1 ( s - t) \eqsim m_2 \pi, \\ 2 m_1 m_2 t \eqsim m_1 m_2 (s + t) - m_2 m_1 ( s - t) \eqsim m_2 \pi. \end{align*} Thus, if |$m_2$| is even, we have |$s = t_{l^{\prime}}^{(\boldsymbol{m})}$| and |$t = t_{l}^{(\boldsymbol{m})}$| for some |$l,l^{\prime} \in{{\mathbb{Z}}}$|⁠. On the other hand, if |$m_2$| is odd, we obtain |$s = t_{l^{\prime}+\frac 12}^{(\boldsymbol{m})}$| and |$t = t_{l+\frac 12}^{(\boldsymbol{m})}$|⁠. If |$m_2$| is even, we can conclude the following: if |$t \in [0, 2\pi )$| and |$t \neq t_{l}^{(\boldsymbol{m})}$| for some |$l \in \{0, \ldots , 2 m_1 m_2 -1\}$|⁠, then |$t$| is the only element of |$[0,2 \pi )$| in |$\mathcal{A}^{(\boldsymbol{m})}(t)$| and the statement (iii) of the proposition is proven. If |$t \in [0, 2\pi )$| and |$t = t_{l}^{(\boldsymbol{m})}$|⁠, |$l\in \{0,\ldots ,2m_1m_2-1\}$| and |$l \not \equiv 0 \mod m_1$|⁠, then |$s \in [0,2\pi )$| given by (2.4), with |$v = -1$| satisfies both conditions in (2.3) for |$v = -1$|⁠. Further, because of the second condition in (2.3), |$s \not \eqsim t$|⁠. Note that the particular choice of the numbers |$a$| and |$b$| from Bézout’s lemma does not influence equation (2.4) so that if |$\boldsymbol{m}$| is fixed and |$v = -1$|⁠, then |$s$| is uniquely determined by |$t$|⁠. Since |$t=t_{l}^{(\boldsymbol{m})}$|⁠, the so constructed |$s$| can also be written as |$s = t_{l^{\prime}}^{(\boldsymbol{m})}$|⁠, with some |$l^{\prime}\in \{0,\ldots ,2m_1m_2-1\}$|⁠, |$l^{\prime} \not \equiv 0 \mod m_1$| and |$l^{\prime} \neq l$|⁠. In total, we can conclude that |$\mathcal{A}^{(\boldsymbol{m})}(t) = \{s,t\}$| and, thus, the statement (ii) of the proposition. If |$m_2$| is odd, we obtain statements (ii)′ and (iii)′ in an analogous way. From the findings in Proposition 2.1 we see that an even frequency number |$m_2$| leads to a slightly different setup of intersection points than an odd |$m_2$|⁠. In this article, we will focus on the case that |$m_2$| is an even number. In this case the nodes \begin{equation} {\textbf{LS}}^{(\boldsymbol{m})}_{\alpha} = \left\{\,\boldsymbol{\ell}^{(\boldsymbol{m})}_{\alpha}\big(t^{(\boldsymbol{m})}_{l}\big)\,|\, l\in \big\{0,\ldots,2m_1m_2-1\big\}\right\} \end{equation} (2.5) contain the two poles of the sphere and give a simple characterization of all self-intersection points of the Lissajous curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$|⁠. Corollary 2.3 Let |$\textrm{gcd}(\boldsymbol{m}) = 1$| and |$m_2$| even. Then, |$\textbf{LS}^{(\boldsymbol{m})}_{\alpha }$| is the set of all self-intersection points of the closed curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t)$|⁠, |$t \in [0,2\pi )$|⁠. |$\textbf{LS}^{(\boldsymbol{m})}_{\alpha }$| contains |$m_2 (m_1 -1) + 2$| points on the sphere |${\mathbb{S}}^2$|⁠, including both poles that are traversed |$m_2$| times, and |$m_2(m_1-1)$| nonpolar double points that are traversed |$2$| times by the curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }(t)$| as |$t$| varies from |$0$| to |$2\pi $|⁠. 3. Characterizing spherical Lissajous nodes In addition to the description given in Corollary 2.3, we can characterize the intersection points of the Lissajous curves also as the union of two interlacing rectangular grids in spherical coordinates. The construction for this second characterization can be performed for general frequencies |$\boldsymbol{m} = (m_1, m_2) \in{{\mathbb{N}}}^{2}$|⁠, where |$m_2$| is even. If |$m_1$| and |$m_2$| are not relatively prime the so obtained nodes can also be interpreted in terms of Lissajous curves. This relation will be discussed at the end of this section. To describe the spherical Lissajous nodes we introduce the index set \begin{equation} \textbf{I}^{(\boldsymbol{m})} = \left\{ \,(i_1, i_2) \in{{\mathbb{N}}}_0^{2}\ \left|\begin{array}{ll} & 0\leq i_{1}\leq m_{1}, \; 0 \leq i_2 < 2 m_2, \\ & i_2 < m_2 \ \textrm{if}\ i_{1}\ \in \{0,m_1\}, \\ & i_{1} + i_{2}\ \textrm{is even} \end{array}\right.\!\!\! \right\}. \end{equation} (3.1) The set |$\textbf{I}^{(\boldsymbol{m})}$| is a disjoint union |$\textbf{I}^{(\boldsymbol{m})} = \textbf{I}^{(\boldsymbol{m})}_0 \cup \textbf{I}^{(\boldsymbol{m})}_1$| of the two sets \begin{equation} \textbf{I}^{(\boldsymbol{m})}_0 = \big\{\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})} \ | \ i_1,\ i_2\ \textrm{are even} \, \big\}, \quad \textbf{I}^{(\boldsymbol{m})}_1 = \big\{\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})} \ | \ i_1,\ i_2\ \textrm{are odd} \, \big\}. \end{equation} (3.2) For |$\boldsymbol{i} = (i_1, i_2) \in \textbf{I}^{(\boldsymbol{m})}$| we obtain a relation to spherical coordinates by introducing the latitudinal and longitudinal angles \begin{equation*} \theta^{(m_1)}_{i_1} = \frac{i_1}{m_1} \pi \in [0,\pi], \qquad \varphi^{(m_2)}_{i_2} = \frac{i_2}{m_2} \pi \in [0,2\pi).\end{equation*} The set of nodes on the sphere |${\mathbb{S}}^2$| corresponding to these spherical coordinates is given by \begin{equation} \textbf{LS}^{(\boldsymbol{m})} =\left\{\, \boldsymbol{x}^{(\boldsymbol{m})}_{\boldsymbol{i}}\,\left|\,\boldsymbol{i}\in \textbf{I}^{(\boldsymbol{m})} \right.\right\}, \end{equation} (3.3) with the points |$\boldsymbol{x}^{(\boldsymbol{m})}_{\boldsymbol{i}} \in{\mathbb{S}}^2$| defined by \begin{equation*}\boldsymbol{x}^{(\boldsymbol{m})}_{\boldsymbol{i}} = \left(\sin \big(\theta^{(m_1)}_{i_1}\big) \cos \big(\varphi^{(m_2)}_{i_2}\big), \sin \big(\theta^{(m_1)}_{i_1}\big) \sin \big(\varphi^{(m_2)}_{i_2}\big), \cos \big(\theta^{(m_1)}_{i_1}\big)\right).\end{equation*} The cardinality of the set |$\textbf{I}^{(\boldsymbol{m})}$| in (3.1) can be determined from the simple structure of the sets |$\textbf{I}^{(\boldsymbol{m})}_0$| and |$\textbf{I}^{(\boldsymbol{m})}_1$| in (3.2) (see also Fig. 2). We have \begin{equation*} \#\textbf{I}^{(\boldsymbol{m})}_{0} = m_1 m_2 / 2, \quad \#\textbf{I}^{(\boldsymbol{m})}_{1} = m_1 m_2 / 2, \end{equation*} Fig. 2. View largeDownload slide Illustration of the index sets |$\textbf{I}^{(\boldsymbol{m})}$|⁠, |$\textbf{I}^{(\boldsymbol{m})}_0$| and |$\textbf{I}^{(\boldsymbol{m})}_1$|⁠, as well as |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$| as defined in (3.1), (3.2) and (3.4), respectively. The black marks describe |$\textbf{I}^{(\boldsymbol{m})}_0$|⁠, the white marks |$\textbf{I}^{(\boldsymbol{m})}_1$| and the union gives |$\textbf{I}^{(\boldsymbol{m})}$|⁠. According to (3.3), the indices |$\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}$| with |$i_1 = 0$| and |$i_1 = m_1$| describe the north and south pole of |${\mathbb{S}}^2$|⁠, respectively. For a one-to-one relation at the poles we consider the reduced subset |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})} \subset \textbf{I}^{(\boldsymbol{m})}$|⁠. The square marks on the left and right boundary describe the indices not contained in |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$|⁠. Fig. 2. View largeDownload slide Illustration of the index sets |$\textbf{I}^{(\boldsymbol{m})}$|⁠, |$\textbf{I}^{(\boldsymbol{m})}_0$| and |$\textbf{I}^{(\boldsymbol{m})}_1$|⁠, as well as |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$| as defined in (3.1), (3.2) and (3.4), respectively. The black marks describe |$\textbf{I}^{(\boldsymbol{m})}_0$|⁠, the white marks |$\textbf{I}^{(\boldsymbol{m})}_1$| and the union gives |$\textbf{I}^{(\boldsymbol{m})}$|⁠. According to (3.3), the indices |$\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}$| with |$i_1 = 0$| and |$i_1 = m_1$| describe the north and south pole of |${\mathbb{S}}^2$|⁠, respectively. For a one-to-one relation at the poles we consider the reduced subset |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})} \subset \textbf{I}^{(\boldsymbol{m})}$|⁠. The square marks on the left and right boundary describe the indices not contained in |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$|⁠. and, thus, \begin{equation*}\#\textbf{I}^{(\boldsymbol{m})} =\#\textbf{I}^{(\boldsymbol{m})}_0 +\#\textbf{I}^{(\boldsymbol{m})}_1 = m_1 m_2. \end{equation*} All points |$\boldsymbol{x}^{(\boldsymbol{m})}_{\boldsymbol{i}}$| with |$i_1 = 0$| describe the north pole of |${\mathbb{S}}^2$| and all |$\boldsymbol{x}^{(\boldsymbol{m})}_{\boldsymbol{i}}$| with |$i_1 = m_1$| the south pole. Therefore, the cardinality of |$\textbf{LS}^{(\boldsymbol{m})}$| is smaller than |$\#\textbf{I}^{(\boldsymbol{m})}$|⁠. A simple counting gives |$\#\textbf{LS}^{(\boldsymbol{m})} = (m_1 -1) m_2 +2$|⁠. In order to have a one-to-one correspondence between indices and elements in |$\textbf{LS}^{(\boldsymbol{m})}$| we introduce the following subsets of |$\textbf{I}^{(\boldsymbol{m})}$|⁠: \begin{equation} \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})} = \big\{\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})} \ | \ i_2 \in \{0,1\} \ \textrm{if} \ i_1 \in \{0, m_1\}\big\}. \end{equation} (3.4) Clearly |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})} \subset \textbf{I}^{(\boldsymbol{m})} $| and |$\# \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})} = \# \textbf{LS}^{(\boldsymbol{m})} = (m_1 -1) m_2 +2$|⁠. Moreover, every node in |$\textbf{LS}^{(\boldsymbol{m})}$| can now be described in a unique way by an index |$\boldsymbol{i} \in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$|⁠. In particular, we have \begin{equation*}{\textbf{LS}}^{(\boldsymbol{m})} =\big\{\, \boldsymbol{x}^{(\boldsymbol{m})}_{\boldsymbol{i}}\,\big|\,\boldsymbol{i}\in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}\big\}.\end{equation*} As a basis for the interpolation on the sphere, we will use a double Fourier basis that is not continuous at the poles of the sphere. Therefore, it makes sense to formulate the interpolation theory first in terms of the larger index set |$\textbf{I}^{(\boldsymbol{m})}$|⁠. In the second step, we will then reduce the problem to the subset |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$| and the corresponding Lissajous node points |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. The reason for the halved number of elements at the left and right boundary in |$\textbf{I}^{(\boldsymbol{m})}$| is a glide reflection symmetry of the used Fourier basis. This symmetry will play an important role when we discuss the implementation of the interpolation scheme. The following more technical result provides an identification of the index set |$\textbf{I}^{(\boldsymbol{m})}$| with a class decomposition of the product set |$H^{(\boldsymbol{m})}\times{R}^{(\boldsymbol{m})}$|⁠, where the sets |$H^{(\boldsymbol{m})}$| and |${R}^{(\boldsymbol{m})}$| are given as \begin{equation*} H^{(\boldsymbol{m})}=\{0,\ldots, 2 m_{1} m_{2}/g -1\} \qquad \textrm{and}\qquad R^{(\boldsymbol{m})}= \{0,\ldots,g-1\}.\end{equation*} Here, |$g = \textrm{gcd}(\boldsymbol{m})$| denotes again the greatest common divisor of the integers |$m_1$| and |$m_2$|⁠. This result will provide us with the link between the nodes |$\textbf{LS}^{(\boldsymbol{m})}$| and the involved generating Lissajous curves. Proposition 3.1 Let |$\boldsymbol{m} \in{{\mathbb{N}}}^2$| and |$m_2$| be even. (a) For all |$(l,\rho )\in H^{(\boldsymbol{m})}\times R^{(\boldsymbol{m})}$|⁠, there exists an |$\boldsymbol{i}\in \textbf{I}^{(\boldsymbol{m})}$| and a |$v \in \{-1,1\}$| such that \begin{align} i_{1} &\equiv v l \mod 2 m_1,\qquad\qquad\qquad\qquad \end{align} (3.5) \begin{align} i_{2} &\equiv l - 2\rho \frac{m_2}{g} - \frac{1-v}{2} m_2 \mod 2 m_2. \end{align} (3.6) The number |$v \in \{-1,1\}$| and the element |$\boldsymbol{i}\in \textbf{I}^{(\boldsymbol{m})}$| are uniquely determined by (3.5) and (3.6). In this way, a function |$\boldsymbol{i}^{(\boldsymbol{m})}:\,H^{(\boldsymbol{m})}\times R^{(\boldsymbol{m})}\to \textbf{I}^{(\boldsymbol{m})}$| is well defined by \begin{equation*} \boldsymbol{i}^{(\boldsymbol{m})}(l,\rho)=\boldsymbol{i}. \end{equation*} (b) |$\boldsymbol{i}^{(\boldsymbol{m})}(l,\rho )\in \textbf{I}^{(\boldsymbol{m})}_{0}$| if and only if |$l$| is even, and |$\boldsymbol{i}^{(\boldsymbol{m})}(l,\rho )\in \textbf{I}^{(\boldsymbol{m})}_{1}$| if and only if |$l$| is odd. (c) For all |$\boldsymbol{i}\in \textbf{I}^{(\boldsymbol{m})}$| we have |$\#\{(l,\rho )\in H^{(\boldsymbol{m})} \times R^{(\boldsymbol{m})}\,|\,\boldsymbol{i}^{(\boldsymbol{m})}(l,\rho )=\boldsymbol{i}\,\} = 2$|⁠. Proof. We start with statement (a). For |$l \in H^{(\boldsymbol{m})}$| we can find an integer |$0\leq i_{1}\leq m_{1}$| and a |$v \in \{-1,1\}$| such that (3.5) is satisfied. Clearly, the number |$i_1$| is uniquely determined by this condition, whereas |$v \in \{-1,1\}$| is only uniquely determined if |$l \not \equiv 0 \mod 2 m_1$| and |$l \not \equiv m_1 \mod 2 m_1$|⁠. Further, for |$(l,\rho )$| and |$v \in \{-1,1\}$| given by (3.5), equation (3.6) gives a uniquely determined integer |$0\leq i_{2}< 2m_{2}$| in the case that |$l$| is not divisible by |$m_1$| or |$2 m_1$|⁠. In the case that |$l \equiv 0 \mod 2 m_1$| or |$l \equiv m_1 \mod 2 m_1$|⁠, condition (3.6) provides a unique integer |$0\leq i_{2}< m_{2}$| by determining at the same time the value of |$v \in \{-1,1\}$|⁠. Since |$m_2$| is even, we have |$i_{1} \equiv l\equiv i_{2} \mod 2$|⁠. This implies |$\boldsymbol{i}\in \textbf{I}^{(\boldsymbol{m})}$| and, therefore, statement (a). Statement (b) follows also directly from (3.5), (3.6) and the definition in (3.2). We finally turn to statement (c). If |$\boldsymbol{i}\in \textbf{I}^{(\boldsymbol{m})}$| and |$v\in \{-1,1\}$|⁠, then the integers |$a_{1}=v i_1$| and |$a_2 = i_2 + m_2 (1-v)/2$| are uniquely determined by |$\boldsymbol{i}\in \textbf{I}^{(\boldsymbol{m})}$|⁠, |$v\in \{-1,1\}$| and satisfy |$a_1 \equiv a_2 \mod 2$|⁠. Since |$m_1$| and |$m_2/g$| are relatively prime the Chinese remainder theorem yields a unique number |$l\in \{0,\ldots ,2m_1m_2/g-1\}$| such that \begin{equation*} v i_{1} \equiv l \mod 2m_1, \qquad i_{2} + m_2 (1-v)/2 \equiv l \mod 2m_2/g. \end{equation*} Now, we can find also a uniquely determined |$\rho \in{R}^{(\boldsymbol{m})}$| such that (3.6) holds. Thus, since both choices of |$v$| give distinct elements |$(l,\rho )$|⁠, statement (c) is shown. A simple consequence of this proposition is the following description of the nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. Corollary 3.2 Let |$\boldsymbol{m} \in{{\mathbb{N}}}^2$| and |$m_2$| even. Then \begin{equation*}{\textbf{LS}}^{(\boldsymbol{m})} = \bigcup_{\rho = 0}^{g-1} \left\{\,\boldsymbol{\ell}^{(\boldsymbol{m})}_{2 \rho/m_2} \big(t^{(\boldsymbol{m})}_{l}\big)\ | \ l\in \{0,\ldots,2m_1m_2/g-1\} \,\right\},\end{equation*} where |$t^{(\boldsymbol{m})}_{l}$| are the equidistant sampling points given in (1) and |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{\alpha }$| the Lissajous curves introduced in (1.1). In particular, if |$m_1$| and |$m_2$| are relatively prime, then |$\textbf{LS}^{(\boldsymbol{m})}$| corresponds to the self-intersection points |$\textbf{LS}^{(\boldsymbol{m})}_0$| of the curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_{0}$| given in (2.5). Proof. By definition of the Lissajous curve and the sampling points |$t^{(\boldsymbol{m})}_{l}$|⁠, we have \begin{equation*} \textstyle \boldsymbol{\ell}^{(\boldsymbol{m})}_{2 \rho/m_2} (t^{(\boldsymbol{m})}_{ l}) = \left( \sin \left( \frac{ l \pi}{m_1} \right) \cos \left( \frac{l \pi}{m_2} - \frac{2 \rho}{m_2} \pi \right), \ \sin \left(\frac{ l \pi}{m_1} \right) \sin \left(\frac{l \pi}{m_2} - \frac{2 \rho}{m_2} \pi \right), \ \cos \left( \frac{ l \pi}{m_1}\right) \right).\end{equation*} Now applying Proposition 3.1 we can find |$\boldsymbol{i}\in \textbf{I}^{(\boldsymbol{m})}$| and |$v \in \{-1,1\}$| such that the relations (3.5) and (3.6) are satisfied. In particular, this implies \begin{align*} \boldsymbol{\ell}^{(\boldsymbol{m})}_{2 \rho /m_2} (t^{(\boldsymbol{m})}_{ l}) &= \textstyle \left( \sin \left( \frac{ v i_1 \pi}{m_1} \right) \cos \left( \frac{ i_2 \pi}{m_2} -\frac{1-v}{2} \pi \right), \ \sin \left( \frac{ v i_1 \pi}{m_1} \right) \sin \left(\frac{ i_2 \pi}{m_2} -\frac{1-v}{2} \pi \right), \ \cos \left( \frac{ v i_1 \pi}{m_1}\right) \right) \\ &= \textstyle \left( \sin \left( \frac{ i_1 \pi}{m_1} \right) \cos \left( \frac{ i_2 \pi}{m_2} \right), \ \sin \left( \frac{ i_1 \pi}{m_1} \right) \sin \left(\frac{ i_2 \pi}{m_2} \pi \right), \ \cos \left( \frac{ i_1 \pi}{m_1}\right) \right) = \boldsymbol{x}_{\boldsymbol{i}}^{(\boldsymbol{m})}. \end{align*} Going these steps back, we get the reverse implication: if |$\boldsymbol{x}_{\boldsymbol{i}}^{(\boldsymbol{m})}$| is given, we can fix |$v \in \{-1,1\}$| and obtain by Proposition 3.1 a unique pair |$(l,\rho )$| such that |$\boldsymbol{x}_{\boldsymbol{i}}^{(\boldsymbol{m})} = \boldsymbol{\ell }^{(\boldsymbol{m})}_{2 \rho /m_2} (t^{(\boldsymbol{m})}_{ l})$|⁠. If |$m_1$| and the even |$m_2$| are relatively prime, Corollary 3.2 provides the second attempt to characterize the self-intersection points of the spherical Lissajous curves mentioned at the beginning of this section. If |$m_1$| and |$m_2$| are not relatively prime, it states that |$\textbf{LS}^{(\boldsymbol{m})}$| can be generated by time equidistant samples of at most |$g$| different Lissajous curves. Two examples of node sets |$\textbf{LS}^{(\boldsymbol{m})}$| in which |$m_1$| and |$m_2$| are not relatively prime are illustrated in Fig. 3. Fig. 3. View largeDownload slide Illustration of the nodes |$\boldsymbol{\textrm{LS}}^{(\boldsymbol{m})}$| and the union of their generating curves for |$m_1 = m_2$|⁠. Fig. 3. View largeDownload slide Illustration of the nodes |$\boldsymbol{\textrm{LS}}^{(\boldsymbol{m})}$| and the union of their generating curves for |$m_1 = m_2$|⁠. 4. Discrete orthogonality structure on |$\textbf{I}^{(\boldsymbol{m})}$| We denote by |$\mathcal{L}(\textbf{I}^{(\boldsymbol{m})})$| the space of all discrete functions on |$\textbf{I}^{(\boldsymbol{m})}$|⁠. For |$\boldsymbol{\gamma }\in{{\mathbb{Z}}}^{2}$|⁠, we consider the functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}\in \mathcal{L}(\textbf{I}^{(\boldsymbol{m})})$| given by \begin{equation} \chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{i}) = \left\{ \begin{array}{ll} \cos(\gamma_{1} i_1 \pi/m_{1})\, {e}^{\imath \gamma_{2} i_2 \pi/m_{2}} & \textrm{if}\ \gamma_2\ \textrm{is even}, \\ \imath\ \sin(\gamma_{1} i_1 \pi/m_{1})\, {e}^{\imath \gamma_{2} i_2 \pi/m_{2}} & \textrm{if}\ \gamma_2\ \textrm{is odd}. \end{array} \right. \end{equation} (4.1) The functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$| are a discretization of the parity-modified Fourier basis that we will discuss in the next section. The goal of this section is to establish a discrete orthogonality of the functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$| on |$\textbf{I}^{(\boldsymbol{m})}$| similar to the discrete orthogonality structure developed for the Lissajous–Chebyshev points in Dencker & Erb (2017a,b). This will be our main technical prerequisite for the proofs of the upcoming interpolation results. We denote the normalized uniform discrete measure on the power set of |$\textbf{I}^{(\boldsymbol{m})}$| by |$\omega ^{(\boldsymbol{m})}$|⁠. It is determined by |$\omega ^{(\boldsymbol{m})}(\{\boldsymbol{i}\})= 1/(m_1m_2)$|⁠. The vector space |$\mathcal{L}(\textbf{I}^{(\boldsymbol{m})})$| endowed with the inner product \begin{equation*} \langle\, f,g \rangle_{\omega^{(\boldsymbol{m})}} = \int f \, \overline{g} \, \mathrm{d}\omega^{(\boldsymbol{m})} = \frac1{m_1m_2}\sum_{\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}} f(\boldsymbol{i}) \overline{g(\boldsymbol{i})}\end{equation*} is a Hilbert space. The corresponding norm is denoted by |$\|\cdot \|_{\omega ^{(\boldsymbol{m})}}$|⁠. Proposition 4.1 Let |$\boldsymbol{m} \in{{\mathbb{N}}}^2$|⁠, |$m_2$| be even and |$\boldsymbol{\gamma }\in{{\mathbb{Z}}}^{2}$|⁠. If |$\int \chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }} \mathrm{d}\omega ^{(\boldsymbol{m})}\neq 0$|⁠, then \begin{equation} \textrm{there exist}\ (h_1,h_2)\ \in\ \mathbb{Z}_0^{2}\ \textrm{with}\ \gamma_{1}=h_{1}m_{1},\ \gamma_{2}=h_{2}m_{2},\ \textrm{and}\ h_1\ +\ h_2\ \equiv\ 0\ \mod\ 2. \end{equation} (4.2) If (4.2) is satisfied, then |$\int \chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }} \mathrm{d}\omega ^{(\boldsymbol{m})}=1$|⁠. In the proof, we use for |$N\in{{\mathbb{N}}}_0$| the well-known trigonometric identity \begin{equation} \sum_{l=0}^N \,{e}^{\imath l \vartheta} = \left\{ \begin{array}{ll} \frac{{e}^{\imath (N+1) \vartheta}-1}{{e}^{\imath \vartheta}-1} \quad & \vartheta\notin 2\pi\mathbb{Z},\\ N+1 & \vartheta \in 2\pi\mathbb{Z}. \end{array} \right. \end{equation} (4.3) Proof. We start with the case |$\gamma _2 \equiv 0 \mod 2$|⁠. Then, using Proposition 3.1, we obtain \begin{align*} \int\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\mathrm{d}\omega^{(\boldsymbol{m})} &= \frac{1}{m_1m_2}\sum_{\boldsymbol{i}\in\textbf{I}^{(\boldsymbol{m})}} \cos(\gamma_{1} i_1 \pi/m_{1}) \,{e}^{\imath \gamma_{2} i_2 \pi/m_{2}} \\ &= \frac{1}{2m_1m_2} \sum_{v \in \{-1,1\}}\sum_{\boldsymbol{i}\in\textbf{I}^{(\boldsymbol{m})}} {e}^{\imath \,( \gamma_{1} v i_1 \pi/m_{1} + \gamma_{2} i_2 \pi/m_{2})}\\ &= \frac{1}{2m_1m_2} \sum_{l \in H^{(\boldsymbol{m})}} \sum_{\rho \in R^{(\boldsymbol{m})}} {e}^{\imath \,( \gamma_{1} l \pi/m_{1} + \gamma_{2} l \pi/m_{2} + 2 \gamma_2 \rho \pi / m_2)}. \end{align*} In view of (4.3), this integral is only different from zero if |$\gamma _{1} /m_{1} + \gamma _{2} /m_{2} \in 2 {{\mathbb{Z}}}$| and |$\gamma _2 /m_2 \in{{\mathbb{Z}}}$| are satisfied. Thus, if we assume that the integral |$\int \chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}\mathrm{d}\omega ^{(\boldsymbol{m})} \neq 0$| then |$\gamma _2 = h_2 m_2$| with some integer |$h_2 \in{{\mathbb{Z}}}$| and |$\gamma _{1} /m_{1} + \gamma _2 / m_2 \in 2 {{\mathbb{Z}}}$|⁠. Thus, also |$\gamma _1$| is of the form |$\gamma _1 = h_1 m_1$| with an integer |$h_1 \in{{\mathbb{Z}}}$| and we further have |$h_1 + h_2 \in 2 {{\mathbb{Z}}}$|⁠. This proves (4.2). On the other hand if (4.2) is satisfied then (4.3) gives \begin{align*} \int\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\mathrm{d}\omega^{(\boldsymbol{m})} &= \frac{1}{2m_1m_2} \sum_{l \in H^{(\boldsymbol{m})}} {e}^{\imath \,l ( \gamma_{1} \pi/m_{1} + \gamma_{2} \pi/m_{2})} \sum_{\rho \in R^{(\boldsymbol{m})}} {e}^{\imath \,( 2 \gamma_2 \pi / g) \rho} = 1. \end{align*} In the case |$\gamma _2 \equiv 1 \mod 2$|⁠, we obtain with Proposition 3.1 \begin{align*} \int\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\mathrm{d}\omega^{(\boldsymbol{m})} &= \frac{1}{m_1m_2} \sum_{\boldsymbol{i}\in\textbf{I}^{(\boldsymbol{m})}} \imath \sin(\gamma_{1} i_1 \pi/m_{1}) \,{e}^{\imath \gamma_{2} i_2 \pi/m_{2}} \\ &= \frac{1}{2m_1m_2} \sum_{v \in \{-1,1\}}\sum_{\boldsymbol{i}\in\textbf{I}^{(\boldsymbol{m})}} {e}^{\imath \,( \gamma_{1} v i_1 \pi/m_{1} + \gamma_{2} i_2 \pi/m_{2} + (1-v)\pi /2 )}\\ &= \frac{1}{2m_1m_2} \sum_{l \in H^{(\boldsymbol{m})}} \sum_{\rho \in R^{(\boldsymbol{m})}} {e}^{\imath \,( \gamma_{1} l \pi/m_{1} + \gamma_{2} l \pi/m_{2} + 2 \gamma_2 \rho \pi /m_2)}. \end{align*} Now, with the same argumentation as above, the fact that this integral does not vanish implies the conditions |$\gamma _1 = h_1 m_1$|⁠, |$\gamma _2 = h_2 m_2$| for some integers |$h_1$| and |$h_2$|⁠, as well as |$h_1 + h_2 \in 2 {{\mathbb{Z}}}$|⁠. Furthermore, if (4.2) is satisfied also in this case we get |$ \int \chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}\mathrm{d}\omega ^{(\boldsymbol{m})} = 1$|⁠. Using the discrete functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$| and Proposition 4.1, we are now going to construct two orthogonal basis systems in the space |$\mathcal{L}(\textbf{I}^{(\boldsymbol{m})})$|⁠. For this, we introduce the spectral index set \begin{equation} \overline{\boldsymbol{\varGamma}}^{(\boldsymbol{m})} = \left\{\,\boldsymbol{\gamma}\in{{\mathbb{N}}} \times{{\mathbb{Z}}} \left|\begin{array}{ll} 1 \leq \gamma_{1} \leq m_{1},\\ \gamma_{1}/m_{1}+ |\gamma_{2}|/m_{2} \leq 1 \end{array}\right. \right\} \cup \left\{\,(0,\gamma_2) \left|\begin{array}{ll} \gamma_2 \in 2 {{\mathbb{Z}}}, \\ |\gamma_{2}| < m_{2} \end{array}\right. \right\}. \end{equation} (4.4) For odd |$\gamma _2$|⁠, we have |$\chi ^{(\boldsymbol{m})}_{(0,\gamma _2)}(\boldsymbol{i}) = 0$| for all |$\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}$|⁠. Hence, these indices are excluded in (4.4). The index set |$\overline{\boldsymbol{\varGamma }}^{(\boldsymbol{m})}$| is in general still too large for our purpose. Some of the functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$|⁠, |$\boldsymbol{\gamma } \in \overline{\boldsymbol{\varGamma }}^{(\boldsymbol{m})}$|⁠, are linearly dependent in |$\mathcal{L}(\textbf{I}^{(\boldsymbol{m})})$|⁠. This linear dependence in |$\overline{\boldsymbol{\varGamma }}^{(\boldsymbol{m})}$| is related to the two sets \begin{align*} \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{U}} &= \left\{ \boldsymbol{\gamma}\in \overline{\boldsymbol{\varGamma}}^{(\boldsymbol{m})} \ | \ \gamma_{1}/m_{1} + \gamma_{2}/m_{2} = 1, \ \gamma_2 \neq 0 \ \right\}, \\ \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}} &= \left\{ \boldsymbol{\gamma}\in \overline{\boldsymbol{\varGamma}}^{(\boldsymbol{m})} \ | \ \gamma_{1}/m_{1} - \gamma_{2}/m_{2} = 1, \ \gamma_2 \neq 0 \ \right\}. \end{align*} In particular, the set given by \begin{equation} \boldsymbol{\varGamma}^{(\boldsymbol{m})} = \overline{\boldsymbol{\varGamma}}^{(\boldsymbol{m})} \setminus \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{U}} \end{equation} (4.5) will soon turn out to be the right spectral index set for our considerations. Two examples of such spectral index sets are given in Fig. 4. Note that the choice of |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{U}}$| over |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{D}}$| in (4.5) is arbitrary and can be also switched for the subsequent results. Also note that if |$m_1$| and |$m_2$| are relatively prime then |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{U}} = \boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{D}} = \emptyset $| and |$\boldsymbol{\varGamma }^{(\boldsymbol{m})} = \overline{\boldsymbol{\varGamma }}^{(\boldsymbol{m})}$|⁠. A simple counting argument gives the following complexities: \begin{equation*} \begin{array}{rlllll} \#\overline{\boldsymbol{\varGamma}}^{(\boldsymbol{m})} &=& m_1 m_2 + g - 1, \qquad\qquad & \#\boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{U}} &=& g-1, \\ \#\boldsymbol{\varGamma}^{(\boldsymbol{m})} &=& m_1 m_2, & \#\boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}} &=& g-1. \end{array} \end{equation*} For the proof of the subsequent theorem, we will use the following product formulas: \begin{align} \chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}} \overline{\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}^{\prime}}} &= \frac{1}{2} \left(\chi^{(\boldsymbol{m})}_{(\gamma_1-\gamma_1^{\prime},\gamma_2 - \gamma_2^{\prime})} + (-1)^{\gamma_2^{\prime}} \chi^{(\boldsymbol{m})}_{{(\gamma_1+\gamma_1^{\prime},\gamma_2 - \gamma_2^{\prime})}} \right), \end{align} (4.6) \begin{align} \chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}} \chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}^{\prime}} &= \frac{1}{2} \left(\chi^{(\boldsymbol{m})}_{(\gamma_1+\gamma_1^{\prime},\gamma_2 + \gamma_2^{\prime})} + (-1)^{\gamma_2^{\prime}} \chi^{(\boldsymbol{m})}_{{(\gamma_1-\gamma_1^{\prime},\gamma_2 + \gamma_2^{\prime})}} \right), \end{align} (4.7) \begin{align} \overline{\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}} &= (-1)^{\gamma_2} \chi^{(\boldsymbol{m})}_{(\gamma_1,-\gamma_2)}. \end{align} (4.8) These identities are satisfied for all |$\boldsymbol{\gamma },\boldsymbol{\gamma }^{\prime} \in{{\mathbb{Z}}}^2$| and can be derived directly from the definition (4.1) of the functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$| using standard identities for the products of two trigonometric functions. Fig. 4. View largeDownload slide Illustration of the index sets |$\boldsymbol{\varGamma }^{(\boldsymbol{m})}$| (black dots) and |$\overline{\boldsymbol{\varGamma }}^{(\boldsymbol{m})}$| (black and white dots). The white dots are the elements of |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{U}}$|⁠. The circled dots indicate the basis functions in (4.9) with doubled norm. Fig. 4. View largeDownload slide Illustration of the index sets |$\boldsymbol{\varGamma }^{(\boldsymbol{m})}$| (black dots) and |$\overline{\boldsymbol{\varGamma }}^{(\boldsymbol{m})}$| (black and white dots). The white dots are the elements of |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{U}}$|⁠. The circled dots indicate the basis functions in (4.9) with doubled norm. Theorem 4.2 Let |$\boldsymbol{m} \in{{\mathbb{N}}}^2$|⁠, |$m_2$| be even. The functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$|⁠, |$\boldsymbol{\gamma }\in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$|⁠, form an orthogonal basis of the |$m_1m_2$|-dimensional inner product space |$(\mathcal{L}(\textbf{I}^{(\boldsymbol{m})}),\langle \,\cdot ,\cdot \,\rangle _{\omega ^{(\boldsymbol{m})}})$|⁠. The norms of the basis functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$| are given as \begin{equation} \|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\|_{\omega^{(\boldsymbol{m})}}^2 = \left\{ \begin{array}{rl} 1,\; & \textrm{if}\quad \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}, \ \gamma_1 \in \{0,m_1\},\\ \frac12,\; & \textrm{if}\quad \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}, \ \gamma_1 \notin \{0,m_1\}. \end{array} \right. \end{equation} (4.9) Proof. We will continuously use the product formula (4.6) in this proof. Therefore, we denote the index vectors on the right-hand side of (4.6) by \begin{equation*} \boldsymbol{\gamma}^+ = (\gamma_1+\gamma_1^{\prime},\gamma_2-\gamma_2^{\prime}) \quad \textrm{and} \quad \boldsymbol{\gamma}^- = (\gamma_1-\gamma_1^{\prime},\gamma_2-\gamma_2^{\prime}).\end{equation*} We assume first that |$\boldsymbol{\gamma },\boldsymbol{\gamma }^{\prime}\in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$| and |$\boldsymbol{\gamma } \neq \boldsymbol{\gamma }^{\prime}$|⁠. We differentiate between two subcases. Case 1:|$\boldsymbol{\gamma } = (m_1,0)$| or |$\boldsymbol{\gamma ^{\prime}} = (m_1,0)$|⁠. Without restriction, we assume that |$\boldsymbol{\gamma } = (m_1,0)$|⁠. Then, we have \begin{equation*} m_1 \leq m_1 + \gamma_1^{\prime} < 2 m_1, \quad 0 < m_1 - \gamma_1^{\prime} \leq m_1, \quad m_2 < \gamma_2^{\prime} < m_2.\end{equation*} This implies that both |$\boldsymbol{\gamma }^+$| and |$\boldsymbol{\gamma }^-$| don’t satisfy the condition (4.2) and therefore, by (4.6), we obtain |$\int \chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }} \overline{\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }^{\prime}}} \mathrm{d}\omega ^{(\boldsymbol{m})} = 0$|⁠. Case 2:|$\boldsymbol{\gamma } \neq (m_1,0)$| and |$\boldsymbol{\gamma ^{\prime}} \neq (m_1,0)$|⁠. Then, based on our assumptions on |$\boldsymbol{\gamma }$| and |$\boldsymbol{\gamma }^{\prime}$|⁠, we obtain the inequalities \begin{equation*} 0 \leq \gamma_1 + \gamma_1^{\prime} < 2 m_1, \quad -m_1 < \gamma_1 - \gamma_1^{\prime} < m_1, \quad 2 m_2 < \gamma_2-\gamma_2^{\prime} < 2 m_2.\end{equation*} For |$\boldsymbol{\gamma }^-$|⁠, the condition (4.2) can only be satisfied if |$\boldsymbol{\gamma }=\boldsymbol{\gamma }^{\prime}$|⁠, which is excluded by the given assumptions. For |$\boldsymbol{\gamma }^+$|⁠, the condition (4.2) is satisfied if |$\boldsymbol{\gamma }=\boldsymbol{\gamma }^{\prime}=(0,0)$| or if |$\gamma _1 + \gamma _1^{\prime} = m_1$| and |$|\gamma _2-\gamma _2^{\prime}| = m_2$| holds true. The first instance can be excluded by the assumption |$\boldsymbol{\gamma } \neq \boldsymbol{\gamma }^{\prime}$|⁠. Also, the second instance can be excluded, since by |$\boldsymbol{\varGamma }^{(\boldsymbol{m})} \subset \overline{\boldsymbol{\varGamma }}^{(\boldsymbol{m})}$| and |$\boldsymbol{\varGamma }^{(\boldsymbol{m})} \cap \boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{U}} = \emptyset $| we have \begin{equation*} \frac{\gamma_1}{m_1} + \frac{\gamma_2}{m_2} + \frac{\gamma_1^{\prime}}{m_1} - \frac{\gamma_2^{\prime}}{m_2} < 2, \quad \frac{\gamma_1}{m_1} - \frac{\gamma_2}{m_2} + \frac{\gamma_1}{m_1} + \frac{\gamma_2}{m_2} < 2. \end{equation*} Thus, by Proposition 4.1 and (4.6) we obtain also for the second case |$\int \chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }} \overline{\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }^{\prime}}} \mathrm{d}\omega ^{(\boldsymbol{m})} = 0$|⁠. In total, we can conclude that the functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$|⁠, |$\boldsymbol{\gamma }\in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$|⁠, are pairwise orthogonal with respect to the inner product |$\langle \,\cdot ,\cdot \,\rangle _{\omega ^{(\boldsymbol{m})}}$|⁠. We now have a look at the norms of the functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$| and consider the case |$\boldsymbol{\gamma } = \boldsymbol{\gamma }^{\prime}$| in (4.6). We get |$\boldsymbol{\gamma }^+ = (2 \gamma _1,0)$| and |$\boldsymbol{\gamma }^- = (0,0)$|⁠. Since |$\boldsymbol{\gamma }\in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$|⁠, we have |$0 \leq 2 \gamma _1 \leq 2 m$|⁠. Therefore, condition (4.2) is always satisfied for |$\boldsymbol{\gamma }^-$| and satisfied for |$\boldsymbol{\gamma }^+$| precisely if |$\gamma _1 \in \{0,m_1\}$|⁠. Based on this observation Proposition 4.1 implies (4.9). Finally, since the functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$|⁠, |$\boldsymbol{\gamma } \in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$|⁠, form an orthogonal system consisting of |$m_1 m_2$| elements and |${\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$| is a space of the same dimension, the functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$|⁠, |$\boldsymbol{\gamma } \in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$|⁠, are an orthogonal basis of |${\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$|⁠. Remark 4.3 Theorem 4.2 holds also true if we replace |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{U}}$| in the definition (4.5) of |$\boldsymbol{\varGamma }^{(\boldsymbol{m})}$| with the counterpart |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{D}}$|⁠. The respective proof is identical. For practical issues it is convenient to have also a real basis for the vector space |${\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$|⁠. For this, we introduce a second set of basis functions as \begin{equation} \chi^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma}}} = \! \left\{ \! \begin{array}{ll} \operatorname{Re} \chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}} & \textrm{if} \quad \begin{array}{l} \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \setminus \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_2 \leq 0, \quad \textrm{or}\\ \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_1 \leq m_1/2, \end{array} \\[5mm] \operatorname{Im} \chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}} &\textrm{if} \quad \begin{array}{l} \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \setminus \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_2> 0, \quad \textrm{or}\\ \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_1 > m_1/2. \end{array} \end{array} \right. \end{equation} (4.10) Theorem 4.4 Let |$\boldsymbol{m} \in{{\mathbb{N}}}^2$|⁠, |$m_2$| be even. The functions |$\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}}$|⁠, |$\boldsymbol{\gamma }\in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$|⁠, are a real orthogonal basis of the inner product space |$({\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})}),\langle \,\cdot ,\cdot \,\rangle _{\omega ^{(\boldsymbol{m})}})$|⁠. The norms for |$\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}}$| read as \begin{equation} \big\|\chi^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma}}}\big\|_{\omega^{(\boldsymbol{m})}}^2 = \left\{ \begin{array}{rl} 1\; &\ \textrm{if} \quad \boldsymbol{\gamma} \in \{(0,0), (m_{1},0) \}, \\[2pt] \frac{1}{2}\; &\ \textrm{if}\ \quad \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \setminus \{(0,0),(m_1,0)\}, \;\ \gamma_1\ =\ 0\ \textrm{or}\ \gamma_2\ =\ 0,\\[2pt] \frac{1}{2}\; &\ \textrm{if}\ \quad \boldsymbol{\gamma} = (m_1/2,-m_2/2), \\[2pt] \frac{1}{4}\; &\ \textrm{for all other}\ \boldsymbol{\gamma}\ \in\ \boldsymbol{\varGamma}^{(\boldsymbol{m})}. \end{array} \right. \end{equation} (4.11) Proof. The functions |$\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}}$| are clearly all real and given as |$\operatorname{Re} \chi ^{(\boldsymbol{m})}_{(\boldsymbol{\gamma })} = \frac 12 (\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }} + \overline{\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}})$| and |$\operatorname{Im} \chi ^{(\boldsymbol{m})}_{(\boldsymbol{\gamma })} = \frac 1{2 \imath } (\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }} - \overline{\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}})$|⁠. Now, based on the formulas (4.6), (4.7) and (4.8) as well as Proposition 4.1, the statements about the orthogonality and the norms of the basis functions can be derived similarly as in Theorem 4.2. As a template for the entire procedure, we calculate the norm |$\|\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}}\|_{\omega ^{(\boldsymbol{m})}}^2$| for the first case given in (4.10), i.e., if |$\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}} = \operatorname{Re} \chi ^{(\boldsymbol{m})}_{(\boldsymbol{\gamma })}$|⁠. Here, we get \begin{align*} \big\|\chi^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma}}}\big\|_{\omega^{(\boldsymbol{m})}}^2 &= \frac14 \int \left| \chi^{(\boldsymbol{m})}_{(\gamma_1, \gamma_2)} + (-1)^{\gamma_2} \chi^{(\boldsymbol{m})}_{(\gamma_1, -\gamma_2)}\right|^2 \mathrm{d}\omega^{(\boldsymbol{m})} \\ &= \frac{1}{8} \int \left( 2 \chi^{(\boldsymbol{m})}_{(0, 0)} +2 \chi^{(\boldsymbol{m})}_{(2 \gamma_1, 0)} + \chi^{(\boldsymbol{m})}_{(0, 2\gamma_2)} + \chi^{(\boldsymbol{m})}_{(0, - 2\gamma_2)} + \chi^{(\boldsymbol{m})}_{(2 \gamma_1, 2\gamma_2)} + \chi^{(\boldsymbol{m})}_{(2 \gamma_1, - 2\gamma_2)} \right) \mathrm{d}\omega^{(\boldsymbol{m})}, \end{align*} where we used the product formulas (4.6), (4.7) and (4.8) to manipulate the function terms in the integral. Next, we check in which cases the condition (4.2) given in Proposition 4.1 is satisfied and determine in this way the value of the norm. If |$\boldsymbol{\gamma } \in \{(0,0), (m_{1},0) \}$|⁠, then condition (4.2) is satisfied for all six spectral functions in the integral and we therefore obtain |$\|\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}}\|_{\omega ^{(\boldsymbol{m})}}^2 = 1$|⁠. If |$\boldsymbol{\gamma } \in \boldsymbol{\varGamma }^{(\boldsymbol{m})} \setminus \{(0,0),(m_1,0)\}$| and |$\gamma _1 = 0$| or |$\gamma _2 = 0$|⁠, then condition (4.2) is satisfied only for three of the given basis functions and we obtain |$\|\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}}\|_{\omega ^{(\boldsymbol{m})}}^2 = \frac 12$|⁠. The same holds true if |$\boldsymbol{\gamma } = (m_1/2,-m_2/2)$|⁠. In the remaining case |$\gamma _1>0$|⁠, |$\gamma _2> 0$|⁠, the condition (4.2) is only satisfied for |$\chi ^{(\boldsymbol{m})}_{(0, 0)}$|⁠, and we thus obtain |$\|\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}}\|_{\omega ^{(\boldsymbol{m})}}^2 = \frac 14$|⁠. 5. Interpolation on spherical Lissajous points We are now ready to set up an interpolation scheme for the Lissajous nodes |$\textbf{LS}^{(\boldsymbol{m})}$| on the sphere |${\mathbb{S}}^2$|⁠. We consider general frequencies |$\boldsymbol{m} \in{{\mathbb{N}}}^2$|⁠, where |$m_2$| is even. For simplicity, we will formulate the interpolation problem in the domain |$[0,\pi ] \times [0,2\pi )$| of the spherical coordinates |$(\theta ,\varphi )$|⁠. By (3.3), the corresponding nodes in spherical coordinates are given as |$(\theta ^{(m_1)}_{i_1}, \varphi ^{(m_2)}_{i_2})$|⁠, |$\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}$|⁠. In case we need a one-to-one correspondence for the poles of |${\mathbb{S}}^2$|⁠, we will restrict ourselves to the index set |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})} \subset \textbf{I}^{(\boldsymbol{m})}$| defined in (3.4). For |$\boldsymbol{\gamma } \in{{\mathbb{Z}}}^2$|⁠, we introduce now the following basis functions in spherical coordinates |$(\theta ,\varphi ) \in [0,\pi ] \times [0,2\pi )$|⁠: \begin{equation*} {X}_{\boldsymbol{\gamma}}(\theta,\varphi) = \left\{ \begin{array}{ll} \cos(\gamma_{1} \theta)\, {e}^{\imath \gamma_{2} \varphi} & \textrm{if}\ \; \gamma_2\ \textrm{is even}, \\ \imath \sin(\gamma_{1} \theta)\,{e}^{\imath \gamma_{2} \varphi} & \textrm{if}\ \gamma_2\ \textrm{is odd}. \end{array} \right. \end{equation*} This Fourier type basis for functions on the unit sphere is exactly the basis introduced in Merilees (1973) and Orszag (1974), and mentioned in the introduction. In the literature Boyd (2000), it is referred to as parity-modified Fourier basis. By |$\varPi $|⁠, we denote the space spanned by all linear combinations of the functions |${X}_{\boldsymbol{\gamma }}$|⁠, |$\boldsymbol{\gamma } \in{{\mathbb{Z}}}^2$|⁠. The interpolation problem we want to solve can be stated as follows: for given data values |$f \in{\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$| we want to find a function |$ P^{(\boldsymbol{m})}_f \in \varPi $| such that \begin{equation} P^{(\boldsymbol{m})}_f \big(\theta^{(m_1)}_{i_1}, \varphi^{(m_2)}_{i_2}\big) = f({\boldsymbol{i}}) \quad \textrm{for all}\quad \boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}. \end{equation} (5.1) In order that (5.1) is uniquely solvable, we have to specify an appropriate subspace of |$\varPi $| for the interpolant |$P^{(\boldsymbol{m})}_f$|⁠. For this, the relation \begin{equation} {X}_{\boldsymbol{\gamma}}\big(\theta^{(m_1)}_{i_1}, \varphi^{(m_2)}_{i_2}\big) = \chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{i}), \qquad \boldsymbol{\gamma}\in{{\mathbb{Z}}}^{2}, \quad \boldsymbol{i}\in \textbf{I}^{(\boldsymbol{m})}, \end{equation} (5.2) between the double Fourier basis |${X}_{\boldsymbol{\gamma }}$| and the discrete orthogonal basis |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$| for |${\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$| plays a crucial role. This relation (5.2) and the results of the previous section motivate the introduction of the interpolation space \begin{equation*}\varPi^{(\boldsymbol{m})} = \textrm{span} \big\{\, {X}_{\boldsymbol{\gamma}} \,\big|\, \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \big\}.\end{equation*} Example 5.1 For even |$m \in{{\mathbb{N}}}$|⁠, consider the frequencies |$\boldsymbol{m} = (m-1,m)$|⁠. Then, |$\varPi ^{(\boldsymbol{m})}$| is exactly the space of all parity-modified basis functions |${X}_{\boldsymbol{\gamma }}$| of total degree |$|\gamma _1| + |\gamma _2| \leq m-1$|⁠. Thus, in this case the points |$\textbf{LS}^{(\boldsymbol{m})}$| can be considered as a spherical analog of the Padua points studied in Bos et al. (2006) and Caliari et al. (2005). For |$\boldsymbol{m} = (m,m)$|⁠, they are a spherical version of the bivariate Morrow–Patterson–Xu points introduced and studied in Morrow & Patterson (1978) and Xu (1996). For general |$\boldsymbol{m} \in{{\mathbb{N}}}^2$|⁠, |$m_2$| even, the theory presented in this paper is a spherical analog of the bivariate interpolation theory based on the nodes of two-dimensional Lissajous curves studied in Dencker & Erb (2017a), Dencker & Erb (2017b), Erb (2016) and Erb et al. (2016). In contrast to the actual work, in the literature usually a tensor-product grid in spherical coordinates is used to construct a spectral interpolation scheme on |${\mathbb{S}}^2$| based on the parity-modified double Fourier basis |${X}_{\boldsymbol{\gamma }}$|⁠, see Boyd (1978), Ganesh et al. (1998), Orszag (1974) and Townsend et al. (2017). The corresponding interpolation spaces are defined as |$\operatorname{span} \{ {X}_{\boldsymbol{\gamma }} \ | \ 0 \leq \gamma _1 \leq m_1, \ |\gamma _2| \leq m_2\}$| by using a rectangular spectral index set. Respective variants are also established for bivariate polynomial interpolation and are sometimes referred to as maximal degree spaces. A comparison between different bivariate interpolation spaces related to total degree and maximal degree spaces can be found in the treatise (Trefethen, 2017). In order to have a one-to-one correspondence between data values on |$\textbf{LS}^{(\boldsymbol{m})}$| and |$\textbf{I}^{(\boldsymbol{m})}$|⁠, we additionally consider the subspaces \begin{align} {\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})}) &= \,f \in \mathcal{L}(\textbf{I}^{(\boldsymbol{m})}) \ | \ f(\boldsymbol{i}) \equiv f(\,\boldsymbol{j}) \quad \textrm{if} \; i_1 = j_1 \in \{0,m_1\}, \\ \varPi^{(\boldsymbol{m})}_{\mathrm{S}} &= P \in \varPi^{(\boldsymbol{m})} \ | \ P(\theta^{(m_1)}_{i_1}, \varphi^{(m_2)}_{i_2}) \equiv P(\theta^{(m_1)}_{j_1}, \varphi^{(m_2)}_{j_2}) \quad \textrm{if} \; i_1 = j_1 \in \{0,m_1\}.\notag \end{align} (5.3) Clearly |${\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})}) \subset{\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$| and |$\dim{\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})}) = \# \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})} = \# \textbf{LS}^{(\boldsymbol{m})} = \dim \varPi ^{(\boldsymbol{m})}_{\mathrm{S}}$|⁠. The data functions |$f \in{\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})})$| are constant at the coordinates |$i_1 = 0$| and |$i_1 = m_1$| corresponding to the poles of the sphere. The space |${\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})})$| can therefore be used to describe all possible data sets on the Lissajous nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. The subspace |$\varPi ^{(\boldsymbol{m})}_{\mathrm{S}} \subset \varPi ^{(\boldsymbol{m})}$| gives all elements |$P \in \varPi ^{(\boldsymbol{m})}$| such that the data set |$p(\boldsymbol{i}) = P(\theta ^{(m_1)}_{i_1}, \varphi ^{(m_2)}_{i_2})$|⁠, |$\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}$|⁠, is contained in |${\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})})$|⁠. Note that, although |$P \in \varPi ^{(\boldsymbol{m})}_{\mathrm{S}}$| satisfies this discrete pole condition, the function |$P \in \varPi ^{(\boldsymbol{m})}_{\mathrm{S}}$| is in general not constant on the entire lines |$\theta = 0$| and |$\theta = \pi $| describing the poles. As a fundamental basis for the interpolation problem (5.1), we introduce for |$\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}$| the Lagrange functions \begin{equation} L^{(\boldsymbol{m})}_{\boldsymbol{i}}(\theta,\varphi) = \frac{1}{m_1m_2} \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} \frac{\overline{\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{i})}}{\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\|_{\omega^{(\boldsymbol{m})}}^2} {X}_{\boldsymbol{\gamma}}(\theta,\varphi), \qquad (\theta,\varphi) \in [0,\pi] \times [0,2\pi). \end{equation} (5.4) For the subset |$\textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$| defined in (3.4) we use the related variant \begin{equation} L^{(\boldsymbol{m})}_{\mathrm{S},\boldsymbol{i}} = \left\{ \begin{array}{ll} L^{(\boldsymbol{m})}_{\boldsymbol{i}} & \textrm{if}\ \boldsymbol{i}\ \in \ \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})},\ i_1\ \neq\ \{0,m_1\}, \\ \sum_{\boldsymbol{j} \in \textbf{I}^{(\boldsymbol{m})}: j_1 = i_1} L^{(\boldsymbol{m})}_{\boldsymbol{j}} & \textrm{if}\ \boldsymbol{i} \in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})},\ i_1\ \in\ \{0,m_1\}. \end{array} \right. \end{equation} (5.5) We can now state our main result. Theorem 5.2 Let |$\boldsymbol{m} \in{{\mathbb{N}}}^2$|⁠, |$m_2$| be even and |$f\in{\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$|⁠. The interpolation problem (5.1) has a unique solution in the polynomial space |$\varPi ^{(\boldsymbol{m})}$| given by \begin{equation*} P^{(\boldsymbol{m})}_{f}(\theta,\varphi) = \sum_{\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}} f({\boldsymbol{i}}) \, L^{(\boldsymbol{m})}_{\boldsymbol{i}}(\theta,\varphi). \end{equation*} The Lagrange functions |$L^{(\boldsymbol{m})}_{\boldsymbol{i}}$|⁠, |$\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}$|⁠, form a basis of the vector space |$\varPi ^{(\boldsymbol{m})}$|⁠. For |$f\in{\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})})$|⁠, the interpolation problem (5.1) has a solution of the form \begin{equation*} P^{(\boldsymbol{m})}_{f}(\theta,\varphi) = \sum_{\boldsymbol{i} \in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}} f({\boldsymbol{i}}) \, L^{(\boldsymbol{m})}_{\mathrm{S},\boldsymbol{i}}(\theta,\varphi). \end{equation*} This solution is unique in the subspace |$\varPi ^{(\boldsymbol{m})}_{\mathrm{S}} \subset \varPi ^{(\boldsymbol{m})}$| spanned by the functions |$L^{(\boldsymbol{m})}_{\mathrm{S},\boldsymbol{i}}$|⁠, |$\boldsymbol{i} \in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$|⁠. Proof. For |$\,\boldsymbol{j} \in \textbf{I}^{(\boldsymbol{m})}$|⁠, let |$\delta _{\boldsymbol{j}}(\boldsymbol{i}) = \delta _{\boldsymbol{i}\boldsymbol{j}}$| be the Dirac function on |$\textbf{I}^{(\boldsymbol{m})}$|⁠. We consider the system |$\delta _{\boldsymbol{j}}$|⁠, |$\boldsymbol{j} \in \textbf{I}^{(\boldsymbol{m})}$|⁠, as an orthogonal basis of the space |${\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$|⁠. By Theorem 4.2, |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$|⁠, |$\boldsymbol{\gamma } \in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$|⁠, is a second orthogonal basis of |${\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$| and we can expand the functions |$\delta _{\boldsymbol{j}}$|⁠, |$\boldsymbol{j} \in \textbf{I}^{(\boldsymbol{m})}$|⁠, as \begin{equation*}\delta_{\boldsymbol{j}}(\boldsymbol{i}) = \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} \frac{\big\langle\;\! \delta_{\boldsymbol{j}},\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}} \big\rangle_{\omega^{(\boldsymbol{m})}}}{\big\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\big\|_{\omega^{(\boldsymbol{m})}}^2} \chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{i}) = \frac{1}{m_1 m_2} \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} \frac{\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{i}) \overline{\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{j})}}{\big\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\big\|_{\omega^{(\boldsymbol{m})}}^2}.\end{equation*} Evaluating the Lagrange function |$L^{(\boldsymbol{m})}_{\boldsymbol{j}}$|⁠, |$\boldsymbol{j} \in \textbf{I}^{(\boldsymbol{m})}$|⁠, at the points |$(\theta ,\varphi ) = (\theta ^{(m_1)}_{i_1}, \varphi ^{(m_2)}_{i_2})$|⁠, |$\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}$|⁠, and using the identity (5.2), we obtain \begin{equation*}L^{(\boldsymbol{m})}_{\boldsymbol{j}}\big(\theta^{(m_1)}_{i_1},\varphi^{(m_2)}_{i_2}\big) = \frac{1}{m_1m_2} \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} \frac{\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{i}) \overline{\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{j})}}{\big\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\big\|_{\omega^{(\boldsymbol{m})}}^2}.\end{equation*} Thus, |$L^{(\boldsymbol{m})}_{\boldsymbol{j}}(\theta ^{(m_1)}_{i_1},\varphi ^{(m_2)}_{i_2}) = \delta _{\boldsymbol{j}}(\boldsymbol{i})$| and for |$f \in{\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$| the function |$P^{(\boldsymbol{m})}_{f}$| satisfies the interpolation condition (5.1). Furthermore, the mapping |$f \to P^{(\boldsymbol{m})}_{f}$| is an injective homomorphism from |${\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$| into the space |$\varPi ^{(\boldsymbol{m})}$|⁠. Since the dimension |${\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$| coincides with the dimension of |$\varPi ^{(\boldsymbol{m})}$|⁠, this homomorphism is indeed an automorphism and the functions |$L^{(\boldsymbol{m})}_{\boldsymbol{j}}$|⁠, |$\boldsymbol{j} \in \textbf{I}^{(\boldsymbol{m})}$|⁠, form a basis of |$\varPi ^{(\boldsymbol{m})}$|⁠. Finally, we see that for |$f$| in the subspace |${\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})})$| the corresponding function |$P^{(\boldsymbol{m})}_{f}$| is a linear combination of the Lagrange functions |$L^{(\boldsymbol{m})}_{\mathrm{S},\boldsymbol{j}}$|⁠, |$\boldsymbol{j} \in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$|⁠. Since the dimensions of the two subspaces coincide, we get also uniqueness here. As in the discrete case, we want to establish the same result also for a real-valued basis. To this end we define for |$\boldsymbol{\gamma } \in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$| the functions \begin{equation*} {X}_{{\mathcal{R}},{\boldsymbol{\gamma}}}(\theta,\varphi) = \left\{ \begin{array}{ll} \cos(\gamma_{1} \theta) \cos (|\gamma_{2}| \varphi ) & \textrm{if} \quad \begin{array}{l} \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \setminus \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_2 \leq 0, \gamma_2 \ \textrm{even} \quad \textrm{or}\\ \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_1 \leq m_1/2, \gamma_2 \ \textrm{even}, \end{array} \\[5mm] \sin(\gamma_{1} \theta) \sin (|\gamma_{2}| \varphi ) & \textrm{if} \quad \begin{array}{l} \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \setminus \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_2 \leq 0, \gamma_2 \ \textrm{odd} \quad \textrm{or}\\ \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_1 \leq m_1/2, \gamma_2 \ \textrm{odd}, \end{array} \\[5mm] \cos(\gamma_{1} \theta) \sin (\gamma_{2} \varphi ) & \textrm{if} \quad \begin{array}{l} \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \setminus \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_2> 0, \gamma_2 \ \textrm{even}\quad \textrm{or}\\ \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_1 > m_1/2, \gamma_2 \ \textrm{even}, \end{array} \\[5mm] \sin(\gamma_{1} \theta) \cos (\gamma_{2} \varphi ) & \textrm{if} \quad \begin{array}{l} \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \setminus \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_2 > 0, \gamma_2 \ \textrm{odd}\quad \textrm{or}\\ \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_1 > m_1/2, \gamma_2 \ \textrm{odd}. \end{array} \end{array} \right. \end{equation*} Evaluating the functions |${X}_{{\mathcal{R}},{\boldsymbol{\gamma }}}$| at the spherical coordinates |$(\theta ^{(m_1)}_{i_1}, \varphi ^{(m_2)}_{i_2})$| and comparing it with the definition given in (4.10), we obtain the identity \begin{equation*} {X}_{{\mathcal{R}},{\boldsymbol{\gamma}}}\big(\theta^{(m_1)}_{i_1}, \varphi^{(m_2)}_{i_2}\big) = \chi^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma}}}(\boldsymbol{i}) \qquad \textrm{for}\ \boldsymbol{\gamma}\ \in\ \boldsymbol{\varGamma}^{(\boldsymbol{m})}\ \textrm{and}\,\, \boldsymbol{i}\ \in\ \textbf{I}^{(\boldsymbol{m})}. \end{equation*} Based on our experience with the complex-valued basis, it makes sense to introduce the interpolation spaces as \begin{equation*}\varPi^{(\boldsymbol{m})}_{{\mathcal{R}}} = \textrm{span} \big\{\, {X}_{{\mathcal{R}},{\boldsymbol{\gamma}}} \,\big|\, \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \big\}\end{equation*} and the Lagrange functions |$L^{(\boldsymbol{m})}_{\mathcal{R},\boldsymbol{i}}(\theta ,\varphi )$| as \begin{equation*} L^{(\boldsymbol{m})}_{\mathcal{R},\boldsymbol{i}}(\theta,\varphi) = \frac{1}{m_1m_2} \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} \frac{\chi^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma}}}(\boldsymbol{i})}{\|\chi^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma}}}\|_{\omega^{(\boldsymbol{m})}}^2} {X}_{{\mathcal{R}},{\boldsymbol{\gamma}}}(\theta,\varphi) \quad (\theta,\varphi) \in [0,\pi] \times [0,2\pi). \end{equation*} The corresponding reduced subspace |$\varPi ^{(\boldsymbol{m})}_{\mathcal{R},\mathrm{S}}$| and Lagrange functions |$L^{(\boldsymbol{m})}_{\mathcal{R},\mathrm{S},\boldsymbol{i}}$| are defined in the same way as in (5.3) and (5.5), respectively. In analogy to Theorem 5.2, we get the following result. Theorem 5.3 Let |$\boldsymbol{m} \in{{\mathbb{N}}}^2$|⁠, |$m_2$| be even, and |$f\in{\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$|⁠. The interpolation problem (5.1) has a unique solution in the space |$\varPi ^{(\boldsymbol{m})}_{{\mathcal{R}}}$| given by the function \begin{equation*} P^{(\boldsymbol{m})}_{\mathcal{R},\,f}(\theta,\varphi) = \sum_{\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}} f({\boldsymbol{i}}) \, L^{(\boldsymbol{m})}_{\mathcal{R},\boldsymbol{i}}(\theta,\varphi). \end{equation*} The Lagrange functions |$L^{(\boldsymbol{m})}_{\mathcal{R},\boldsymbol{i}}$|⁠, |$\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}$|⁠, form a basis of the vector space |$\varPi ^{(\boldsymbol{m})}_{{\mathcal{R}}}$|⁠. If |$f\in{\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})})$|⁠, the interpolation problem (5.1) has a solution of the form \begin{equation*} P^{(\boldsymbol{m})}_{\mathcal{R},\,f}(\theta,\varphi) = \sum_{\boldsymbol{i} \in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}} f({\boldsymbol{i}}) \, L^{(\boldsymbol{m})}_{\mathcal{R},\mathrm{S},\boldsymbol{i}}(\theta,\varphi). \end{equation*} This solution is unique in the subspace |$\varPi ^{(\boldsymbol{m})}_{\mathcal{R},\mathrm{S}} \subset \varPi ^{(\boldsymbol{m})}_{{\mathcal{R}}}$| spanned by the functions |$L^{(\boldsymbol{m})}_{\mathcal{R},\mathrm{S},\boldsymbol{i}}$|⁠, |$\boldsymbol{i} \in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$|⁠. Remark 5.4 In the discrete setting both basis systems |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$| and |$\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}}$|⁠, |$\boldsymbol{\gamma } \in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$| span the same space |${\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$|⁠. This is different in the continuous setup. Here, we have |$\varPi ^{(\boldsymbol{m})} = \varPi ^{(\boldsymbol{m})}_{{\mathcal{R}}}$| if and only if |$m_1$| and |$m_2$| are relatively prime. If |$m_1$| and |$m_2$| are not relatively prime, then the real basis functions |${X}_{{\mathcal{R}},{\boldsymbol{\gamma }}}$| for |$\gamma \in \boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{D}}$| are linear combinations of complex basis functions |${X}_{\boldsymbol{\gamma }}$| in which the indices |$\gamma $| are contained in both sets |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{D}}$| and |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{U}}$|⁠. 6. Implementation of the interpolation scheme The interpolating function |$P^{(\boldsymbol{m})}_{f}$| can be computed efficiently by using fast Fourier techniques. To this end, we expand |$P^{(\boldsymbol{m})}_{f}$| in the basis |${X}_{\boldsymbol{\gamma }}$| as \begin{equation} P^{(\boldsymbol{m})}_{f}(\theta,\varphi) = \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} c_{\boldsymbol{\gamma}}(\,f) {X}_{\boldsymbol{\gamma}}(\theta,\varphi).\end{equation} (6.1) In this way, once the coefficients |$c_{\boldsymbol{\gamma }}(\,f)$| are calculated, it only remains to evaluate the sum in (6.1). By Theorem 5.2 and definition (5.4) we have the following decomposition: \begin{equation*}P^{(\boldsymbol{m})}_{f}(\theta,\varphi) = \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} \frac{1}{\big\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\big\|_{\omega^{(\boldsymbol{m})}}^2} \left( \frac{1}{m_1 m_2}\sum_{\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}} f({\boldsymbol{i}}) \overline{\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{i})} \right) {X}_{\boldsymbol{\gamma}}(\theta,\varphi).\end{equation*} Since the functions |${X}_{\boldsymbol{\gamma }}$| form a basis of |$\varPi ^{(\boldsymbol{m})}$|⁠, we immediately obtain the following. Corollary 6.1 For |$\,f\in{\mathcal{L}}(\textbf{I}^{(\boldsymbol{m})})$|⁠, the uniquely determined coefficients |$c_{\boldsymbol{\gamma }}(f)$| in the expansion (6.1) are given by \begin{equation*}c_{\boldsymbol{\gamma}}(\,f)=\textstyle \frac{1}{\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\|_{\omega^{(\boldsymbol{m})}}^2}\,\big\langle\, f,\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}} \big\rangle_{\omega^{(\boldsymbol{m})}}.\end{equation*} Calculation of the coefficients |$c_{\boldsymbol{\gamma }}(\,f)\ \ $| Based on the formula in Corollary 6.1, the coefficients |$c_{\boldsymbol{\gamma }}(\,f)$| can be computed using a two-dimensional Fourier transform on the finite abelian group |${{\mathbb{Z}}} / (2 m_1) \times{{\mathbb{Z}}} / (2 m_2)$|⁠. We identify this group with \begin{equation*} {\boldsymbol{J}}^{(\boldsymbol{m})} = \big\{ \boldsymbol{i} \in{{\mathbb{Z}}}^2 \ | \ 0 \leq i_1 \leq 2 m_1 -1, \ 0 \leq i_2 \leq 2 m_2 -1 \big\}. \end{equation*} We introduce a flip operator on |${\boldsymbol{J}}^{(\boldsymbol{m})}$| by defining |$\boldsymbol{i}^* = (2 m_1 - i_1 \mod 2 m_1, i_2 + m_2 \mod 2 m_2)$| for |$\boldsymbol{i} \in{\boldsymbol{J}}^{(\boldsymbol{m})}$|⁠. Using this operator, we can extend a function |$f$| on |$\textbf{I}^{(\boldsymbol{m})}$| to a function |$g$| on the whole group |${\boldsymbol{J}}^{(\boldsymbol{m})}$| by setting \begin{equation*} g(\boldsymbol{i}) = \frac{1}{2 m_1m_2} \left\{ \begin{array}{rl} f(\boldsymbol{i}), \quad &\quad \textrm{if}\; \boldsymbol{i}\in\textbf{I}^{(\boldsymbol{m})}, \,\\ f(\boldsymbol{i}^*), \quad &\quad \textrm{if}\; \boldsymbol{i}^*\in\textbf{I}^{(\boldsymbol{m})}, \,\\ 0, \quad &\quad \textrm{otherwise}. \end{array} \right. \end{equation*} The computation of the coefficient |$c_{\boldsymbol{\gamma }}(\,f)$| can now be reduced to the calculation of the Fourier transform |$\hat{g}$| of |$g$| on |${\boldsymbol{J}}^{(\boldsymbol{m})}$| by using the identity \begin{align}c_{\boldsymbol{\gamma}}(\,f)&=\textstyle \frac{1}{\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\|_{\omega^{(\boldsymbol{m})}}^2}\,\big\langle\, f,\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}} \big\rangle_{\omega^{(\boldsymbol{m})}} =\textstyle \frac{1}{\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\|_{\omega^{(\boldsymbol{m})}}^2}\,\frac{1}{m_1 m_2}\sum_{\boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}} f(\boldsymbol{i}) \overline{\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{i})} \notag \\ & = \textstyle \frac{1}{\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\|_{\omega^{(\boldsymbol{m})}}^2}\,\sum_{\boldsymbol{i} \in{\boldsymbol{J}}^{(\boldsymbol{m})}} g(\boldsymbol{i}) \,{e}^{-\imath \gamma_{1} i_1 \pi/m_{1}} \,{e}^{-\imath \gamma_{2} i_2 \pi/m_{2}} = \left\{ \begin{array}{rl} 2 \hat{g}(\boldsymbol{\gamma}), & \textrm{if}\, \boldsymbol{\gamma}\in\boldsymbol{\varGamma}^{(\boldsymbol{m})}, \, \gamma_1 \notin \{0,m_1\}, \,\\ \hat{g}(\boldsymbol{\gamma}), & \textrm{if}\, \boldsymbol{\gamma}\in\boldsymbol{\varGamma}^{(\boldsymbol{m})}, \, \gamma_1 \in \{0,m_1\}. \end{array} \right. \end{align} (6.2) The computation of the discrete Fourier transform |$\hat{g}(\boldsymbol{\gamma }) = \sum _{\boldsymbol{i} \in{\boldsymbol{J}}^{(\boldsymbol{m})}} g(\boldsymbol{i}) \,{e}^{-\imath \gamma _{1} i_1 \pi /m_{1}} \, {e}^{-\imath \gamma _{2} i_2 \pi /m_{2}}$| can be executed very efficiently in |${\mathcal{O}}(m_1m_2 \log (m_1m_2))$| arithmetic operations using standard algorithms for the fast Fourier transform. The values for |$\|\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}\|_{\omega ^{(\boldsymbol{m})}}^2$| are taken from (4.9). Remark 6.2 The invariance of the function |$g$| under the flip operator, i.e., |$g(\boldsymbol{i}^*) = g(\boldsymbol{i})$|⁠, implies for the Fourier domain the identity |$\hat{g}(\boldsymbol{\gamma }) = (-1)^{\gamma _2} \hat{g}(2m_1 - \gamma _1, \gamma _2)$| for all |$\boldsymbol{\gamma }$| in the dual group (we identify it here also with |${\boldsymbol{J}}^{(\boldsymbol{m})}$|⁠). This glide reflection symmetry of |$g$| and the flip operator are already used in the first publications studying double Fourier series on the sphere (Merilees, 1973; Boer & Steinberg, 1975). In numerical software packages as for instance in Chebfun (Driscoll et al., 2014), this symmetry is used to obtain sparse tensor-product approximations of functions on the sphere (Townsend et al., 2017). In Townsend et al. (2017), the symmetry of |$g$| is called block-mirror-centrosymmetric (BMC) structure. Calculation of the real coefficients |$c_{\mathcal{R},\boldsymbol{\gamma }}(\,f)\ \ $| Also for the real-valued basis |${X}_{{\mathcal{R}},{\boldsymbol{\gamma }}}$| we get an expansion for the interpolating polynomial |$P^{(\boldsymbol{m})}_{\mathcal{R},\,f}$| of the form \begin{equation*} P^{(\boldsymbol{m})}_{\mathcal{R},\,f}(\theta,\varphi) = \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} c_{\mathcal{R},\boldsymbol{\gamma}}(\,f) {X}_{{\mathcal{R}},{\boldsymbol{\gamma}}}(\theta,\varphi),\end{equation*} in which the expansion coefficients |$c_{\mathcal{R},\boldsymbol{\gamma }}(\,f)$| are given by \begin{equation*}c_{\mathcal{R},\boldsymbol{\gamma}}(\,f)=\textstyle \frac{1}{\big\|\chi^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma}}}\big\|_{\omega^{(\boldsymbol{m})}}^2}\,\big\langle\,f,\chi^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma}}} \big\rangle_{\omega^{(\boldsymbol{m})}}.\end{equation*} The calculation of the expansion coefficients can be conducted efficiently using the formula \begin{equation*} c_{\mathcal{R},\boldsymbol{\gamma}}(\,f) = \frac{1}{\|\chi^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma}}}\|_{\omega^{(\boldsymbol{m})}}^2} \left\{ \begin{array}{ll} \operatorname{Re} \hat{g}(\boldsymbol{\gamma}) & \textrm{if} \quad \begin{array}{l} \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \setminus \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_2 \leq 0, \quad \textrm{or}\\ \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_1 \leq m_1/2, \end{array} \\[5mm] - \operatorname{Im} \hat{g}(\boldsymbol{\gamma}) & \textrm{if} \quad \begin{array}{l} \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})} \setminus \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_2> 0, \quad \textrm{or}\\ \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \gamma_1 > m_1/2. \end{array} \end{array} \right. \end{equation*} This formula can be verified as in (6.2) using the real basis (4.10) instead of the complex functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$|⁠. The values |$\|\chi ^{(\boldsymbol{m})}_{{\mathcal{R}},{\boldsymbol{\gamma }}}\|_{\omega ^{(\boldsymbol{m})}}^2$| are explicitly known from (4.11). Calculation of averaged interpolants Instead of using the expansion (6.1), it is sometimes more convenient to implement the more symmetric expansion \begin{equation*} {P}^{(\boldsymbol{m})}_{\mathcal{A},\,f}(\theta,\varphi) = \sum_{\boldsymbol{\gamma} \in \overline{\boldsymbol{\varGamma}}^{(\boldsymbol{m})}} c_{\mathcal{A},\boldsymbol{\gamma}}(\,f) {X}_{\boldsymbol{\gamma}}(\theta,\varphi),\end{equation*} in which the coefficients |$c_{\mathcal{A},\boldsymbol{\gamma }}(\,f)$| are given by \begin{equation*} c_{\mathcal{A},\boldsymbol{\gamma}}(\,f) = \left\{ \begin{array}{ll} c_{\boldsymbol{\gamma}}(\,f)/2 & \textrm{if} \ \boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{U}} \cup \boldsymbol{\varGamma}^{(\boldsymbol{m}),\mathrm{D}}, \\ c_{\boldsymbol{\gamma}}(\,f) & \textrm{for all other}\ \boldsymbol{\gamma}\ \in\ \overline{\boldsymbol{\varGamma}}^{(\boldsymbol{m})}. \end{array}\right.\end{equation*} In this way, it is not necessary to make a choice between |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{U}}$| and |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{D}}$| in order to define the interpolation space. For |$\gamma \in \boldsymbol{\varGamma }^{(\boldsymbol{m}),\mathrm{U}}$| we have |$\chi ^{\boldsymbol{m}}_{\boldsymbol{\gamma }} = (-1)^{\gamma _2} \chi ^{\boldsymbol{m}}_{(m_1 - \gamma _1,-m_2 + \gamma _2)}$| on |$\textbf{I}^{(\boldsymbol{m})}$|⁠, and therefore also |$c_{\boldsymbol{\gamma }}(\,f) = (-1)^{\gamma _2} c_{(m_1 - \gamma _1, - m_2 + \gamma _2)}(\,f)$|⁠. This guarantees that |${P}^{(\boldsymbol{m})}_{\mathcal{A},\,f}$| is also a solution of the interpolation problem (5.1), although in a different space than |$\varPi ^{(\boldsymbol{m})}$|⁠. A similar strategy is of course also possible for the real-valued basis |${X}_{{\mathcal{R}},{\boldsymbol{\gamma }}}$|⁠. Averaged interpolation spaces of this type were originally used for the Morrow–Patterson–Xu points in Harris (2010) and Xu (1996). A more detailed discussion of this averaging related to multivariate interpolation on Lissajous–Chebyshev nodes can be found in Dencker & Erb (2017a). The inverse transform From the coefficients |$c_{\boldsymbol{\gamma }}(\,f)$|⁠, |$\boldsymbol{\gamma } \in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$|⁠, the values |$f \in \mathcal{L}(\textbf{I}^{(\boldsymbol{m})})$| can be recovered efficiently by a second discrete Fourier transform. We give a short description of this inverse transform. Using the interpolation condition (5.1) and (5.2), we have \begin{equation*} f({\boldsymbol{i}}) = P^{(\boldsymbol{m})}_f \big(\theta^{(m_1)}_{i_1}, \varphi^{(m_2)}_{i_2}\big) = \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} c_{\boldsymbol{\gamma}}(\,f) \chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}(\boldsymbol{i}) \quad \textrm{for all}\quad \boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}. \end{equation*} Defining the discrete function |$h$| on the (dual) group |${\boldsymbol{J}}^{(\boldsymbol{m})}$| as \begin{equation*} h(\boldsymbol{\gamma}) = \left\{ \begin{array}{rl} 2 c_{\boldsymbol{\gamma}}(\,f), \quad &\quad \textrm{if}\; \boldsymbol{\gamma}\in\boldsymbol{\varGamma}^{(\boldsymbol{m})}, \, \gamma_1 \notin \{0,m_1\}, \,\\ (-1)^{\gamma_2} 2 c_{\boldsymbol{\gamma}}(\,f), \quad &\quad \textrm{if}\; (2m_1-\gamma_1,\gamma_2)\in\boldsymbol{\varGamma}^{(\boldsymbol{m})}, \, \gamma_1 \notin \{0,m_1\}, \,\\ c_{\boldsymbol{\gamma}}(\,f), \quad &\quad \textrm{if}\; \boldsymbol{\gamma}\in\boldsymbol{\varGamma}^{(\boldsymbol{m})}, \, \gamma_1 \in \{0,m_1\}, \,\\ 0, \quad &\quad \textrm{otherwise}, \end{array} \right. \end{equation*} we obtain from the relation above and the definition (4.1) of the basis functions |$\chi ^{(\boldsymbol{m})}_{\boldsymbol{\gamma }}$| the following discrete Fourier sum: \begin{equation*} f({\boldsymbol{i}}) = \sum_{\boldsymbol{\gamma} \in{\boldsymbol{J}}^{(\boldsymbol{m})}} h(\boldsymbol{\gamma}) \,{e}^{\imath \gamma_{1} i_1 \pi/m_{1}} \,{e}^{\imath \gamma_{2} i_2 \pi/m_{2}} \quad \textrm{for all}\quad \boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}. \end{equation*} In this way, the function |$f \in \mathcal{L}(\textbf{I}^{(\boldsymbol{m})})$| can be recovered by applying a discrete adjoint Fourier transform to |$h$| on |${\boldsymbol{J}}^{(\boldsymbol{m})}$|⁠. As for the computation of the coefficients |$c_{\boldsymbol{\gamma }}(f)$|⁠, this adjoint transform can be executed efficiently in |${\mathcal{O}}(m_1m_2 \log (m_1m_2))$| arithmetic operations. Note that by (6.2) the function |$h$| corresponds to |$\hat{g}$| on |$\boldsymbol{\varGamma }^{(\boldsymbol{m})}$| and |$\{\boldsymbol{\gamma } \in{\boldsymbol{J}}^{(\boldsymbol{m})} \ | \ (2m_1-\gamma _1,\gamma _2) \in \boldsymbol{\varGamma }^{(\boldsymbol{m})}\}$|⁠, but in general not on the entire set |${\boldsymbol{J}}^{(\boldsymbol{m})}$|⁠. 7. Numerical condition and convergence of the interpolation scheme We provide a mathematical description of central properties of the given interpolation scheme, as its numerical condition number, its convergence rates and the behavior at the poles of |${\mathbb{S}}^2$|⁠. The interpolation spaces |$\varPi ^{(\boldsymbol{m})}$| and |$\varPi ^{(\boldsymbol{m})}_{{\mathcal{R}}}$| are spanned by a double Fourier basis with a glide-reflection symmetry. Our strategy is therefore to use the theory of multivariate Fourier series to derive the pursued properties. We consider interpolating functions |${P^{(\boldsymbol{m})}_{f}}$| in which the data |$f$| is given by the samples of a continuous function on the sphere. In particular, if |${\mathscr{f}}(\theta ,\varphi)$| describes a continuous function on |${\mathbb{S}^2}$| in spherical coordinates, we have \begin{equation} f(\boldsymbol{i}) = \mathscr{f}\big(\theta^{(m_1)}_{i_1}, \varphi^{(m_2)}_{i_2}\big) \quad \textrm{for}\quad \boldsymbol{i} \in \textbf{I}^{(\boldsymbol{m})}. \end{equation} (7.1) Clearly, |$f \in{\mathcal{L}}_{\mathrm{S}}(\textbf{I}^{(\boldsymbol{m})})$| and Theorem 5.2 give a unique interpolant |$P^{(\boldsymbol{m})}_{f}$| in |$\varPi ^{(\boldsymbol{m})}_{\mathrm{S}} \subset \varPi ^{(\boldsymbol{m})}$|⁠. Behavior at the poles of the sphere We can describe a continuous function |$\mathscr{f}$| on |${\mathbb{S}}^2$| as a continuous function in spherical coordinates |$(\theta ,\varphi ) \in [0,\pi ] \times [0,2\pi ]$| by using topological identifications at the boundaries. The corresponding function space is given as \begin{equation*} C({\mathbb{S}}^2) = \left\{\mathscr{f} \in C([0,\pi]\times[0,2\pi]) \ \left| \ \begin{array}{lll} \mathrm{(i)} & \mathscr{f}(\theta,0) = \mathscr{f}(\theta,2\pi), & 0 \leq \theta \leq \pi, \\ \mathrm{(ii)} & \mathscr{f}(0,\varphi_1) = \mathscr{f}(0,\varphi_2), & 0 \leq \varphi_1,\varphi_2 \leq 2 \pi, \\ \mathrm{(iii)} &\mathscr{f}(\pi,\varphi_1) = \mathscr{f}(\pi,\varphi_2), & 0 \leq \varphi_1,\varphi_2 \leq 2 \pi. \end{array} \right. \right\} \end{equation*} The parity-modified basis functions |${X}_{\boldsymbol{\gamma }}$|⁠, |$\boldsymbol{\gamma } \in \boldsymbol{\varGamma }^{(\boldsymbol{m})}$|⁠, are in general not contained in |$C({\mathbb{S}}^2)$|⁠. While |${X}_{\boldsymbol{\gamma }} \in C([0,\pi ]\times [0,2\pi ])$| and the periodicity |$\mathrm{(i)}$| are satisfied, the pole conditions |$\mathrm{(ii)}$| and |$\mathrm{(iii)}$| are only satisfied if |$\gamma _2$| is odd. Also the interpolant |$P^{(\boldsymbol{m})}_{f}$| does in general not satisfy the properties |$\mathrm{(ii)}$| and |$\mathrm{(iii)}$| and is therefore not necessarily continuous at the poles of |${\mathbb{S}}^2$|⁠. The condition |$P^{(\boldsymbol{m})}_{f} \in C({\mathbb{S}}^2)$| can be guaranteed only for particular continuous functions |$\mathscr{f}$|⁠. An important example is the space |$\varPi ^{(\boldsymbol{m})} \cap C({\mathbb{S}}^2)$|⁠. Since |$P^{(\boldsymbol{m})}_{f}$| is a projection into |$\varPi ^{(\boldsymbol{m})}$|⁠, we get for |$\mathscr{f} \in C({\mathbb{S}}^2) \cap \varPi ^{(\boldsymbol{m})}$| the identity |$P^{(\boldsymbol{m})}_{f} = \mathscr{f}$| and, thus, |$P^{(\boldsymbol{m})}_{f} \in C({\mathbb{S}}^2)$|⁠. In the next part we will see that the discontinuities of |$P^{(\boldsymbol{m})}_{f}$| at the poles do not affect the global convergence of the interpolation scheme if the function |$\mathscr{f}$| is sufficiently smooth. This guarantees that the interpolant |$P^{(\boldsymbol{m})}_{f}$| and also its derivatives will approximately satisfy conditions |$\mathrm{(ii)}$| and |$\mathrm{(iii)}$| with high accuracy when frequencies |$m_1$| and |$m_2$| get large. In Boyd (1978), such a property at the poles is called a natural boundary condition. Although such a natural condition is sufficient for a lot of applications, there are cases in which the described singularities at the poles result in problems. This pole problem related to the usage of the parity-modified double Fourier basis as well as possible solutions are discussed in Boyd (1978) and Orszag (1974). The Lebesgue constant The operator norm \begin{equation*}\varLambda^{(\boldsymbol{m})} = \sup_{\|\mathscr{f}\|_{\infty} \leq 1} \|P^{(\boldsymbol{m})}_{f}\|_{\infty}, \quad \textrm{with} \;\| \mathscr{f} \|_{\infty} = \sup_{(\theta,\varphi)} |\mathscr{f}(\theta,\varphi)|, \end{equation*} is usually referred to as Lebesgue constant or as absolute condition number of the interpolation problem (5.1). It is an upper bound for the propagation of the error in the uniform norm when constructing the interpolant |$P^{(\boldsymbol{m})}_{f}$| from a continuous function |$\mathscr{f}$|⁠. Theorem 7.1 The Lebesgue constant |$\varLambda ^{(\boldsymbol{m})}$| is bounded by \begin{equation*} \varLambda^{(\boldsymbol{m})} \leq C_{\varLambda} \ln (m_1+1) \ln (m_2+1)\end{equation*} with a constant |$C_{\varLambda }$| independent of |$\boldsymbol{m}$|⁠. Proof. We use the representations (6.1) and (6.2) to rewrite |$P^{(\boldsymbol{m})}_{f}$| in terms of a trigonometric sum. Using the convention |$\hat{g}(-\gamma _1,\gamma _2) = \hat{g}(2 m_1 - \gamma _1,\gamma _2)$| and the glide-reflection symmetry |$\hat{g}(\boldsymbol{\gamma }) = (-1)^{\gamma _2} \hat{g}(2 m_1 - \gamma _1, \gamma _2)$| of |$g$|⁠, we get \begin{align*} P^{(\boldsymbol{m})}_{f}(\theta,\varphi) &= \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m})}} \frac{\hat{g}(\boldsymbol{\gamma})}{\|\chi^{(\boldsymbol{m})}_{\boldsymbol{\gamma}}\|_{\omega^{(\boldsymbol{m})}}^2} {X}_{\boldsymbol{\gamma}}(\theta,\varphi) = \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),*}} \hat{g}(\boldsymbol{\gamma}) \,{e}^{\imath ( \gamma_{1} \theta + \gamma_2 \varphi)} - \hat{g}(m_1,0) \cos(m_1 \theta), \end{align*} where |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),*} = \{ \boldsymbol{\gamma } \in{{\mathbb{Z}}}^2 \ | \ (|\gamma _1|, \gamma _2) \in \boldsymbol{\varGamma }^{(\boldsymbol{m})} \}$| is the symmetric extension of |$\boldsymbol{\varGamma }^{(\boldsymbol{m})}$| from |${{\mathbb{N}}} \times{{\mathbb{Z}}}$| into |${{\mathbb{Z}}}^2$|⁠. For the operator norm |$\varLambda ^{(\boldsymbol{m})} = \sup _{\|\mathscr{f}\|_{\infty } \leq 1} \|P^{(\boldsymbol{m})}_{f}\|_{\infty }$| we get in this way the estimates \begin{align*} \Lambda^{(\boldsymbol{m})} & \leq \sup_{\|\mathscr{f}\|_{\infty} \leq 1} \sup_{(\theta,\varphi)} \left| \sum_{\boldsymbol{\gamma} \in \boldsymbol{\Gamma}^{(\boldsymbol{m}),*}} \sum_{\boldsymbol{i} \in{\boldsymbol{J}}^{(\boldsymbol{m})}} g(\boldsymbol{i}) \,{e}^{-\imath \gamma_{1} (i_1 \pi/m_{1} - \theta) } \,{e}^{-\imath \gamma_{2} (i_2 \pi/m_{2}-\varphi) }\right| + 1 \\ &\leq \sup_{(\theta,\varphi)} \frac{1}{2 m_1 m_2} \sum_{\boldsymbol{i} \in{\boldsymbol{J}}^{(\boldsymbol{m})}} \left| \sum_{\boldsymbol{\gamma} \in \boldsymbol{\Gamma}^{(\boldsymbol{m}),*}} \,{e}^{-\imath \gamma_{1} (i_1 \pi/m_{1} - \theta) } \,{e}^{-\imath \gamma_{2} (i_2 \pi/m_{2}-\varphi) }\right| + 1 \\ &\leq C \int_0^{2\pi} \!\!\! \int_0^{2\pi} \left| \sum_{\boldsymbol{\gamma} \in \boldsymbol{\Gamma}^{(\boldsymbol{m}),*}} {e}^{\imath (\gamma_{1} \theta + \gamma_{2} \varphi) }\right| \mathrm{d}\theta \, \mathrm{d}\varphi + 1. \end{align*} The last transition from a discrete sum to a continuous double integral with a constant |$C>0$| independent of |$\boldsymbol{m}$| is a twofold application of a Marcinkiewicz–Zygmund inequality, see Zygmund (2002, X, Theorem 7.10). The double integral in the last line is known as Fourier–Lebesgue constant of the set |$\boldsymbol{\varGamma }^{(\boldsymbol{m}),*}$|⁠. Taking apart a missing subset |$\{(0,\gamma _2) \ | \ |\gamma _2| \leq m_2, \ \text{$\gamma_2$ odd}\}$|⁠, the Fourier–Lebesgue constants of such sets were studied in Dencker et al. (2017). From the derivations in Dencker et al. (2017, Section 2) (the Lebesgue constant of the missing set is bounded by |$C \ln (m_2 + 1)$|⁠), we get \begin{equation*} \int_0^{2\pi} \!\!\! \int_0^{2\pi} \left| \sum_{\boldsymbol{\gamma} \in \boldsymbol{\varGamma}^{(\boldsymbol{m}),*}} {e}^{\imath (\gamma_{1} \theta + \gamma_{2} \varphi)} \right| \mathrm{d}\theta \, \mathrm{d}\varphi \leq C^{\prime} \ln (m_1+1) \ln (m_2+1),\end{equation*} and, thus, the statement of the theorem. Uniform convergence of the interpolation scheme For a continuous function |$\mathscr{f}$| in spherical coordinates, we denote by |$P^*$| the best approximation in the space |$\varPi ^{(\boldsymbol{m})}$| given by |$ \|\mathscr{f} - P^*\|_{\infty } = \min _{P \in \varPi ^{(\boldsymbol{m})}} \|\mathscr{f} - P\|_{\infty }$|⁠. Using the fact that the interpolation operator |$\mathscr{f} \to P^{(\boldsymbol{m})}_{f}$| reproduces |$P^* \in \varPi ^{(\boldsymbol{m})}$|⁠, together with the bound in Theorem 7.1, we obtain \begin{align*} \|\mathscr{f} - P^{(\boldsymbol{m})}_{f} \|_{\infty} &\leq \|\mathscr{f} - P^* \|_{\infty} + \| P^* - P^{(\boldsymbol{m})}_{f} \|_{\infty} \\ & \leq (\varLambda^{(\boldsymbol{m})} + 1) \|\,f - P^* \|_{\infty} = (C_{\varLambda}+1) \ln (m_1+1) \ln (m_2+1) \|\mathscr{f} - P^* \|_{\infty}. \end{align*} If |$\mathscr{f}$| is |$s$| times continuously differentiable on the sphere, the best error |$\|\mathscr{f} - P^* \|_{\infty }$| can be estimated using a multivariate version of Jackson’s inequality for trigonometric functions, as for instance described in Timan (1963, Section 5.3). As a consequence, we obtain the error estimate \begin{equation} \|\mathscr{f} - P^{(\boldsymbol{m})}_{f}\|_{\infty} \leq C_{\mathscr{f},s} \ln(m_1+1) \ln(m_2+1) \left( \frac{1}{m_1^s} + \frac{1}{m_2^s}\right),\end{equation} (7.2) with a constant |$C_{\mathscr{f},s}$| that depends on |$\mathscr{f}$| and the smoothness |$s$| but not on |$\boldsymbol{m}$|⁠. This kind of error estimate is typical for a multitude of spectral interpolation and approximation methods, and yields a fast uniform convergence of the interpolant provided the original function |$\mathscr{f}$| is smooth. For multivariate polynomial interpolation on Lissajous nodes in the hypercube |$[-1,1]^{\mathsf{d}}$| similar derivations can, for instance, be found in Dencker et al. (2017)and Erb (2016). For a tensor product spectral collocation scheme on the sphere |${\mathbb{S}}^2$|⁠, a corresponding result is provided in Ganesh et al. (1998). 8. Clenshaw–Curtis quadrature formula Using the expansion (6.1), we can easily derive a Clenshaw–Curtis-type interpolatory quadrature formula on the sphere based on function evaluations on |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. In spherical coordinates the area element on the sphere |${\mathbb{S}}^2$| is given by |$\sin \theta \, \mathrm{d} \theta\, \mathrm{d} \varphi $|⁠. Then, the tensor product structure of the basis functions |${X}_{\boldsymbol{\gamma }}$| yields \begin{equation*} \frac{1}{4\pi} \int_0^{2\pi} \int_0^{\pi} {X}_{\boldsymbol{\gamma}}(\theta,\varphi) \sin \theta \, \mathrm{d} \theta\, \mathrm{d} \varphi = \left\{ \begin{array}{ll} \frac{1}{1-\gamma_1^2} & \textrm{if}\ \gamma_2 = 0, \gamma_1\ \textrm{even}, \\ 0 & \textrm{otherwise}.\end{array} \right.\end{equation*} Here, we used the fact that |$\int _{0}^{2\pi } e^{\imath \gamma _2 \varphi } \mathrm{d} \varphi = 2 \pi \delta _{\gamma _2,0}$| and that |$\int _0^{\pi } \cos ( \gamma _1 \theta ) \sin \theta \,\mathrm{d} \theta = \frac{2}{1+\gamma _1^2}$| if |$\gamma _1$| is even and zero otherwise. For the interpolating function |$P^{(\boldsymbol{m})}_{f}$| with the expansion (6.1) we therefore obtain the formula \begin{equation*} \frac{1}{4\pi} \int_0^{2\pi} \int_0^{\pi} P^{(\boldsymbol{m})}_{f}(\theta,\varphi) \sin \theta \, \mathrm{d} \theta \,\mathrm{d} \varphi = \sum_{k = 0}^{\lfloor m_1 /2 \rfloor} \frac{c_{(2k,0)}(f)}{1-4 k^2}. \end{equation*} The coefficients |$c_{(2k,0)}(\,f)$| on the right-hand side depend only on the values |$f(\boldsymbol{i})$|⁠, |$\boldsymbol{i} \in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}$|⁠, which, by (7.1), are linked to function values at the Lissajous nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. This formula can therefore be considered as a Clenshaw–Curtis-type quadrature rule on |${\mathbb{S}}^2$| at the nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. By construction, it is an exact quadrature rule for all functions in |$\varPi ^{(\boldsymbol{m})}$|⁠. Since |$c_{(k,0)}(\,f) = c_{\mathcal{R},(k,0)}(\,f)$|⁠, the same formula holds also true with |$P^{(\boldsymbol{m})}_{\mathcal{R},\,f}$| as an interpolating function. Fig. 5. View largeDownload slide Rotation estimation on the sphere based on sample measurements on the nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. Fig. 5. View largeDownload slide Rotation estimation on the sphere based on sample measurements on the nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. Table 2 Approximation error |$ \| P_{f}^{(\boldsymbol{m})} - \mathscr{f} \|_{\infty }$| for the smooth function |$\mathscr{f}$| in (9.2) |$\boldsymbol{m}$| |$\# \textbf{LS}^{(\boldsymbol{m})}$| |$ \| P_{f}^{(\boldsymbol{m})} - \mathscr{f} \|_{\infty }$| |$\boldsymbol{m}$| |$\# \textbf{LS}^{(\boldsymbol{m})}$| |$ \| P_{f}^{(\boldsymbol{m})} -\mathscr{f} \|_{\infty }$| (3, 4) 16 0.89150031122784 (23, 24) 576 0.00000145422054 (7, 8) 64 0.17505763622726 (27, 28) 784 0.00000003014093 (11, 12) 144 0.01926746577677 (31, 32) 1024 0.00000000047887 (15, 16) 256 0.00126029913111 (35, 36) 1296 0.00000000000604 (19, 20) 400 0.00005152647682 (39, 40) 1600 0.00000000000006 |$\boldsymbol{m}$| |$\# \textbf{LS}^{(\boldsymbol{m})}$| |$ \| P_{f}^{(\boldsymbol{m})} - \mathscr{f} \|_{\infty }$| |$\boldsymbol{m}$| |$\# \textbf{LS}^{(\boldsymbol{m})}$| |$ \| P_{f}^{(\boldsymbol{m})} -\mathscr{f} \|_{\infty }$| (3, 4) 16 0.89150031122784 (23, 24) 576 0.00000145422054 (7, 8) 64 0.17505763622726 (27, 28) 784 0.00000003014093 (11, 12) 144 0.01926746577677 (31, 32) 1024 0.00000000047887 (15, 16) 256 0.00126029913111 (35, 36) 1296 0.00000000000604 (19, 20) 400 0.00005152647682 (39, 40) 1600 0.00000000000006 View Large Table 2 Approximation error |$ \| P_{f}^{(\boldsymbol{m})} - \mathscr{f} \|_{\infty }$| for the smooth function |$\mathscr{f}$| in (9.2) |$\boldsymbol{m}$| |$\# \textbf{LS}^{(\boldsymbol{m})}$| |$ \| P_{f}^{(\boldsymbol{m})} - \mathscr{f} \|_{\infty }$| |$\boldsymbol{m}$| |$\# \textbf{LS}^{(\boldsymbol{m})}$| |$ \| P_{f}^{(\boldsymbol{m})} -\mathscr{f} \|_{\infty }$| (3, 4) 16 0.89150031122784 (23, 24) 576 0.00000145422054 (7, 8) 64 0.17505763622726 (27, 28) 784 0.00000003014093 (11, 12) 144 0.01926746577677 (31, 32) 1024 0.00000000047887 (15, 16) 256 0.00126029913111 (35, 36) 1296 0.00000000000604 (19, 20) 400 0.00005152647682 (39, 40) 1600 0.00000000000006 |$\boldsymbol{m}$| |$\# \textbf{LS}^{(\boldsymbol{m})}$| |$ \| P_{f}^{(\boldsymbol{m})} - \mathscr{f} \|_{\infty }$| |$\boldsymbol{m}$| |$\# \textbf{LS}^{(\boldsymbol{m})}$| |$ \| P_{f}^{(\boldsymbol{m})} -\mathscr{f} \|_{\infty }$| (3, 4) 16 0.89150031122784 (23, 24) 576 0.00000145422054 (7, 8) 64 0.17505763622726 (27, 28) 784 0.00000003014093 (11, 12) 144 0.01926746577677 (31, 32) 1024 0.00000000047887 (15, 16) 256 0.00126029913111 (35, 36) 1296 0.00000000000604 (19, 20) 400 0.00005152647682 (39, 40) 1600 0.00000000000006 View Large 9. Application to rotation estimation on the sphere As a final application, we show how the interpolation scheme presented in this manuscript can be used to estimate the rotation of a function on the sphere based on sample values at the nodes |$\textbf{LS}^{(\boldsymbol{m})}$|⁠. The algorithm to estimate the Euler angles |$\boldsymbol{\beta }= (\beta _1, \beta _2, \beta _3)$| of the rotation follows the scheme presented in Ullisch (2012, Section 4.1). Let |$\mathscr{f}: {\mathbb{S}}^2 \to{{\mathbb{R}}}$| denote the nonrotated function and |$\mathscr{f}_{\textrm{rot}}(\boldsymbol{x}) = \mathscr{f}(R_{\boldsymbol{\beta }} \boldsymbol{x})$| the rotated function on the sphere, where |$R_{\boldsymbol{\beta }}$| denotes the rotation matrix determined by the three Euler angles |$\boldsymbol{\beta }$|⁠. To estimate |$\boldsymbol{\beta }$|⁠, we consider only the data values |$f(\boldsymbol{i}) = \mathscr{f}(\boldsymbol{x}_{\boldsymbol{i}}^{(\boldsymbol{m})})$| and |$\,f_{\textrm{rot}}(\boldsymbol{i})=\mathscr{f}(R_{\boldsymbol{\beta }} \boldsymbol{x}_{\boldsymbol{i}}^{(\boldsymbol{m})})$| measured along the spherical Lissajous curve |$\boldsymbol{\ell }^{(\boldsymbol{m})}_0$|⁠. As an interpolant for the data |$f$| on the unit sphere, we use the function |$P_{f}^{(\boldsymbol{m})}$| (in Cartesian coordinates). In order to obtain an estimate for the Euler angles |$\boldsymbol{\beta }$| we solve the nonlinear least squares problem \begin{equation} \sum_{\boldsymbol{i} \in \textbf{I}_{\mathrm{S}}^{(\boldsymbol{m})}} \big|\,f_{\textrm{rot}}(\boldsymbol{i}) - P_{f}^{(\boldsymbol{m})}\big(R_{\boldsymbol{\beta}} \, \boldsymbol{x}_{\boldsymbol{i}}^{(\boldsymbol{m})} \big)\big|^2 = \textrm{min}. \end{equation} (9.1) In the example given in Fig. 5 we used a linear combination of two Gaussians \begin{equation} \mathscr{f}(\boldsymbol{x}) = e^{-3(x^2 + y^2 + (z-1)^2)} +e^{-4 ((x-1/\sqrt{2})^2 + (y+1/\sqrt{2})^2+z^2)} \end{equation} (9.2) as a test function. As underlying Lissajous curve we chose |$\boldsymbol{\ell }_{0}^{(15,16)}$|⁠. The nonlinear least squares problem (9.1) was solved iteratively with a damped Gauss–Newton scheme. With initial vector |$\boldsymbol{\beta }_0 = (0,0,0)$|⁠, the solution |$(1.4,0.2,0.9)$| is obtained after |$16$| iterations and with the residual |$2.9 \cdot 10^{-3}$|⁠. This result is illustrated in Fig. 5. The approximation errors for the interpolants are listed in Table 2. Note that in general the functional (9.1) has many local minima and particular care has to be given to the choice of the initial vector. A Matlab code of the presented computational example and the developed interpolation scheme on spherical Lissajous nodes is provided at https://github.com/WolfgangErb/LSphere. Acknowledgements I want to thank both referees very much for their excellent work. Their suggestions helped me a lot to improve and extend this manuscript. References Boer , G. & Steinberg , L. ( 1975 ) Fourier series on spheres . Atmosphere , 13 , 180 – 191 . Google Scholar Crossref Search ADS Bos , L. , Caliari , M. , De Marchi , S. , Vianello , M. & Xu , Y. ( 2006 ) Bivariate Lagrange interpolation at the Padua points: the generating curve approach . J. Approx. Theory , 143 , 15 – 25 . Google Scholar Crossref Search ADS Boyd , J. P. ( 1978 ) The choice of spectral functions on a sphere for boundary and eigenvalue problem: a comparison of Chebyshev, Fourier and associated Legendre expansions . Mon. Wea. Rev. , 106 , 1184 – 1191 . Google Scholar Crossref Search ADS Boyd , J. P. ( 2000 ) Chebyshev and Fourier spectral methods . New York : Dover Publications. Caliari , M. , De Marchi , S. & Vianello , M. ( 2005 ) Bivariate polynomial interpolation on the square at new nodal sets . Appl. Math. Comput. , 165 , 261 – 274 . Costa , A. F. , Yen , Y.-F. & Drangova , M. ( 2010 ) Registering spherical navigators with spherical harmonic expansions to measure three-dimensional rotations in magnetic resonance imaging . Magn. Reson. Imaging, 28 , 185 – 194 . Dencker , P. & Erb , W. ( 2017a ) A unifying theory for multivariate polynomial interpolation on general Lissajous-Chebyshev nodes. arXiv : 1711.00557 ( math.NA ). Dencker , P. & Erb , W. ( 2017b ) Multivariate polynomial interpolation on Lissajous-Chebyshev nodes . J. Appr. Theory , 219 , 15 – 45 . Google Scholar Crossref Search ADS Dencker , P. , Erb , W. , Kolomoitsev , Y. & Lomako , T. ( 2017 ) Lebesgue constants for polyhedral sets and polynomial interpolation on Lissajous-Chebyshev nodes . J. Complexity , 43 , 1 – 27 . Google Scholar Crossref Search ADS Driscoll , T. A. , Hale , N. & Trefethen , L. N. (eds) ( 2014 ) Chebfun Guide . Oxford : Pafnuty Publications . Erb , W. ( 2016 ) Bivariate Lagrange interpolation at the node points of Lissajous curves—the degenerate case . Appl. Math. Comput. , 289 , 409 – 425 . Erb , W. , Kaethner , C. , Ahlborg , M. & Buzug , T. M. ( 2016 ) Bivariate Lagrange interpolation at the node points of non-degenerate Lissajous curves . Numer. Math. , 133 , 685 – 705 . Google Scholar Crossref Search ADS Fornberg , B. ( 1995 ) A pseudospectral approach for polar and spherical geometries . SIAM J. Sci. Comp. , 16 , 1071 – 1081 . Google Scholar Crossref Search ADS Ganesh , M. , Graham , I. & Sivaloganathan , J. ( 1998 ) A new spectral boundary integral collocation method for three-dimensional potential problems . SIAM J. Numer. Anal. , 35 , 778 – 804 . Google Scholar Crossref Search ADS Harris , L. A. ( 2010 ) Bivariate Lagrange interpolation at the Chebyshev nodes . Proc. Am. Math. Soc. , 138 , 4447 – 4453 . Google Scholar Crossref Search ADS Merilees , P. E. ( 1973 ) The pseudospectral approximation applied to the shallow water wave equations on a sphere . Atmosphere , 11 , 13 – 20 . Google Scholar Crossref Search ADS Morrow , C. R. & Patterson , T. N. L. ( 1978 ) Construction of algebraic cubature rules using polynomial ideal theory . SIAM J. Numer. Anal. , 15 , 953 – 976 . Google Scholar Crossref Search ADS Orszag , S. A. ( 1974 ) Fourier series on spheres . Mon. Wea. Rev. , 102 , 56 – 75 . Google Scholar Crossref Search ADS Timan , A. F. ( 1963 ) Theory of Approximation of Functions of a Real Variable, translated by J. Berry . Oxford : Pergamon Press . Townsend , A. , Wilber , H. & Wright , G. ( 2017 ) Computing with functions in spherical and polar coordinates I . The sphere. SIAM J. Sci. Comp. , 38 , C403 – C425 . Google Scholar Crossref Search ADS Trefethen , L. N. ( 2017 ) Multivariate polynomial approximation in the hypercube . Proc. Amer. Math. Soc. , 145 , 4837 – 4844 . Google Scholar Crossref Search ADS Ullisch , M. ( 2012 ) A navigator based rigid body motion correction for magnetic resonance imaging. Dissertation , Technische Hochschule Aachen . Welch , E. B. , Manduca , A. , Grimm , R. C. , Ward , H. A. & Clifford , R. J. J. ( 2002 ) Spherical navigator echoes for full 3D rigid body motion measurement in MRI . Magn. Reson. Med. , 47 , 32 – 41 . Google Scholar Crossref Search ADS PubMed Xu , Y. ( 1996 ) Lagrange interpolation on Chebyshev points of two variables . J. Approx. Theory , 87 , 220 – 238 . Google Scholar Crossref Search ADS Zygmund , A. ( 2002 ) Trigonometric Series , vol. I & II combined, 3rd edn. Cambridge : Cambridge University Press . © The Author(s) 2018. Published by Oxford University Press on behalf of the Institute of Mathematics and its Applications. All rights reserved. This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/open_access/funder_policies/chorus/standard_publication_model)
TI - A spectral interpolation scheme on the unit sphere based on the nodes of spherical Lissajous curves
JF - IMA Journal of Numerical Analysis
DO - 10.1093/imanum/dry090
DA - 2018-12-29
UR - https://www.deepdyve.com/lp/oxford-university-press/a-spectral-interpolation-scheme-on-the-unit-sphere-based-on-the-nodes-D0NLJpFSGx
SP - 1
VL - Advance Article
IS - 
DP - DeepDyve
ER -