---
abstract: 'We discuss a novel sampling theorem on the sphere developed by McEwen & Wiaux recently through an association between the sphere and the torus. To represent a band-limited signal exactly, this new sampling theorem requires less than half the number of samples of other equiangular sampling theorems on the sphere, such as the canonical Driscoll & Healy sampling theorem. A reduction in the number of samples required to represent a band-limited signal on the sphere has important implications for compressive sensing, both in terms of the dimensionality and sparsity of signals. We illustrate the impact of this property with an inpainting problem on the sphere, where we show superior reconstruction performance when adopting the new sampling theorem.'
author:
- |
    Jason D. McEwen, Gilles Puy, Jean-Philippe Thiran, Pierre Vandergheynst,\
    Dimitri Van De Ville and Yves Wiaux Institute of Electrical Engineering, Ecole Polytechnique F[é]{}d[é]{}rale de Lausanne (EPFL),\
    Lausanne 1015, Switzerland;\
    Institute of Bioengineering, Ecole Polytechnique F[é]{}d[é]{}rale de Lausanne (EPFL),\
    Lausanne 1015, Switzerland;\
    Department of Radiology and Medical Informatics, University of Geneva (UniGE),\
    Geneva 1211, Switzerland
bibliography:
- 'bib.bib'
title: Sampling theorems and compressive sensing on the sphere
---

INTRODUCTION
============

The fast Fourier transform[@cooley:1965] (FFT) is arguably the most important and widely used numerical algorithm of our era, rendering the frequency content of signals accessible. Moreover, Shannon’s sampling theorem[@shannon:1949] states that all of the information content of a band-limited continuous signal can be captured through a finite number of samples. Typically, standard Fourier analyses are performed in Euclidean space where Shannon’s theory holds and where FFTs are directly applicable. However, in many applications data are observed on non-Euclidean manifolds, such as the sphere. Fourier analysis is performed on the sphere in the basis of spherical harmonics, which are the eigenfunctions of the spherical Laplacian operator and form the canonical orthonormal basis on the sphere. Sampling theorems and fast algorithms to perform spherical harmonic analyses exist but the field is much less mature that its elder Euclidean sibling.

A novel sampling theorem on the sphere has been developed recently by two of the authors of this article [@mcewen:fssht] (hereafter referred to as the MW sampling theorem). From an information theoretic viewpoint, the fundamental property of any sampling theorem is the number of samples required to capture all of the information content of a band-limited signal. To represent exactly a signal on the sphere band-limited at , all sampling theorems on the sphere require $\order(\elmax^2)$ samples. However, the MW sampling theorem requires $\sim2\elmax^2$ samples only, less than half of the $\sim4\elmax^2$ samples of other equiangular sampling theorems on the sphere, such as the canonical Driscoll & Healy sampling theorem[@driscoll:1994] (hereafter referred to as the DH sampling theorem). Not only is the MW sampling theorem of theoretical interest, particularly regarding the information content of signals on the sphere, but it also has important practical implications in the emerging field of compressive sensing.

The theory of compressive sensing states that it is possible to acquire sparse or compressible signals with fewer samples than standard sampling theorems would suggest [@candes:2006a; @donoho:2006]. In these settings, the ratio of the required number of measurements to the dimensionality of the signal scales linearly with its sparsity [@candes:2006a]. The more efficient sampling of the MW sampling theorem reduces the dimensionality of the signal in the spatial domain, thereby improving the performance of compressive sensing reconstruction on the sphere when compared to alternative sampling theorems.[@mcewen:css2] Furthermore, for sparsity priors defined in the spatial domain, such as signals sparse in the magnitude of their gradient, sparsity is also directly related to the sampling of the signal. For this class of signals, an additional enhancement in compressive sensing reconstruction performance is thus achieved when adopting the MW sampling theorem. [@mcewen:css2]

In this article we first review sampling theorems on the sphere in , focussing on the MW and DH sampling theorems. In we discuss the superior performance achieved when solving compressive sensing problems on the sphere using the MW sampling theorem, as opposed to the DH sampling theorem. We illustrate our arguments with an inpainting problem on the sphere, where we adopt the prior assumption that the signal to be recovered is sparse in the magnitude of its gradient. Simulations are performed, verifying our theoretical arguments. Finally, concluding remarks are made in .

SAMPLING THEOREMS ON THE SPHERE {#sec:sampling_theorems}
===============================

Sampling theorems on the sphere state that all of the information contained in a band-limited signal may be represented by a finite set of samples in the spatial domain. On the sphere, unlike Euclidean space, the number of samples required in the harmonic and spatial domains differ, with different sampling theorems on the sphere requiring a different number of samples in the spatial domain. For an equiangular sampling of the sphere, the DH sampling theorem has become the standard, while the MW sampling theorem has emerged only recently.[^1] The MW sampling theorem achieves a more efficient sampling of the sphere, with a reduction by a factor of two in the number of samples required to represent a band-limited signal.[^2] In this section we outline the harmonic structure of the sphere in the continuous setting, before reviewing concisely the DH and MW sampling theorems.

Harmonic Analysis on the Sphere
-------------------------------

We consider the space of square integrable functions on the sphere $\ltwo(\sphere)$, with the inner product of $f,g\in\ltwo(\sphere)$ defined by $$\innerp{f}{g} = \int_\sphere \dmu{\sas} \: f(\sas) \: g^\cconj(\sas) 
\spcend ,$$ where $\dmu{\sas} = \sin \saa \dx \saa \dx \sab$ is the usual invariant measure on the sphere and $(\sas)$ define spherical coordinates with colatitude $\saa \in [0,\pi]$ and longitude $\sab \in [0,2\pi)$. Complex conjugation is denoted by the superscript ${}^\cconj$.

The spherical harmonic functions form the canonical orthogonal basis for the space of $\ltwo(\sphere)$ functions on the sphere and are defined by $$\shfarg{\el}{\m}{\sas} = 
\sqrt{\frac{2\el+1}{4\pi} \elmfact} \:
\aleg{\el}{\m}{\cos\saa} \:
\exp{\img \m \sab} 
\spcend ,$$ for natural $\el\in\naturals$ and integer $\m\in\integers$, $|\m|\leq\el$, where $\aleg{\el}{\m}{x}$ are the associated Legendre functions. We adopt the Condon-Shortley phase convention, with the $(-1)^\m$ phase factor included in the definition of the associated Legendre functions. The orthogonality and completeness relations for the spherical harmonics read $
\innerp{\shf{\el}{\m}}{\shf{\el\p}{\m\p}}
= 
\kron{\el}{\el\p}
\kron{\m}{\m\p}
$ and $$\sumlm
\shfarg{\el}{\m}{\saa,\sab} \:
\shfargc{\el}{\m}{\saa\p,\sab\p} 
=
\delta(\cos\saa - \cos\saa\p) \:
\delta(\sab - \sab\p)$$ respectively, where $\kron{i}{j}$ is the Kronecker delta symbol and $\delta(x)$ is the Dirac delta function.

Due to the orthogonality and completeness of the spherical harmonics, any square integrable function on the sphere $\f \in \ltwo(\sphere)$ may be represented by its spherical harmonic expansion $$\label{eqn:sht_inverse}
\f(\sas) = 
\sum_{\el=0}^\infty
\sum_{\m=-\el}^\el
\flm \:
\shfarg{\el}{\m}{\sas}
\spcend ,$$ where the spherical harmonic coefficients are given by the usual projection onto each basis function: $$\label{eqn:sht_forward}
\shc{\f}{\el}{\m} =
\innerp{f}{\shf{\el}{\m}}
=
\int_\sphere \dmu{\sas} \: f(\sas) \: \shfargc{\el}{\m}{\sas}
\spcend .$$ Throughout, we consider signals on the sphere band-limited at $\elmax$, that is signals such that $\shc{f}{\el}{\m}=0$, $\forall
\el\geq\elmax$. In this case the summation over  in [Eqn. (\[eqn:sht\_inverse\])]{} may be truncated to $\elmax-1$. We also adopt the convention that harmonic coefficients are defined to be zero for $|\m|>\el$ (which enforces the contraint $|\m| \leq \el$ when summations are interchanged).

The forward and inverse spherical harmonic transforms, given by [Eqn. (\[eqn:sht\_forward\])]{} and [Eqn. (\[eqn:sht\_inverse\])]{} respectively, are exact in the continuous setting. A sampling theorem on the sphere states how to sample a band-limited function $\f(\sas)$ at a finite number of locations, such that all of the information content of the continuous function is captured. Since the frequency domain of a function on the sphere is discrete, the spherical harmonic coefficients describe the continuous function exactly. A sampling theorem thus describes how to exactly recover the spherical harmonic coefficients  of the continuous function from its samples. Consequently, sampling theorems effectively encode (often implicitly) an exact quadrature rule for evaluating the integral of a band-limited function on the sphere.

Driscoll & Healy Sampling Theorem
---------------------------------

The DH sampling theorem[@driscoll:1994] gives an explicit quadrature rule for the evaluation of spherical harmonic coefficients: $$\label{eqn:dh_quad}
\shc{\f}{\el}{\m} =
\sum_{\saai=0}^{2\elmax-1} \:
\sum_{\sabi=0}^{2\elmax-1} \:
\qweightdh(\saaiang) \: \f(\saisang) \: \shfargc{\el}{\m}{\saisang}
\spcend ,$$ where the equiangular sample positions are defined by $\saaiang= \pi
\saai / 2\elmax $, for $\saai = 0, \dotsc, 2\elmax-1$, and $\sabiang =
\pi\sabi / \elmax$, for $\sabi = 0, \dotsc, 2\elmax -1$, giving $\Ndh=(2\elmax-1)2\elmax + 1 \sim 4\elmax^2$ samples on the sphere.[^3] The quadrature used to evaluate [Eqn. (\[eqn:sht\_forward\])]{} is exact for a function band-limited at , with quadrature weights given implicitly by the solution to $$\label{eqn:dh_quad_weight_implicit}
\sum_{\saai=0}^{2\elmax-1} \:
\qweightdh(\saaiang) \: \leg{\el}{\cos \saaiang}
= \frac{2\pi}{\elmax} \: \kron{\el}{0}
\spcend ,
\quad \forall \el < 2\elmax
\spcend .$$ The quadrature weights satisfying [Eqn. (\[eqn:dh\_quad\_weight\_implicit\])]{} are given by $$\qweightdh(\saaiang) =
\frac{2\pi}{\elmax^2} \:
\sin\saaiang \:
\sum_{k=0}^{\elmax-1} \:
\frac{\sin((2k+1)\saaiang)}{2k+1}
\spcend .$$

The exactness of the quadrature rule is proved by considering the sampling distribution of Dirac delta functions defined by $$\label{eqn:dh_sampling_distn}
s(\sas) = 
\sum_{\saai=0}^{2\elmax-1} \:
\sum_{\sabi=0}^{2\elmax-1} \:
\qweightdh(\saaiang) \: 
\delta(\cos\saa - \cos\saaiang) \: \delta(\sab - \sabiang)
\spcend .$$ For quadrature weights satisfying [Eqn. (\[eqn:dh\_quad\_weight\_implicit\])]{}, it can be shown that $\shc{s}{0}{0}=\sqrt{4\pi}$ and $\shc{s}{\el}{\m}=0$ for $0<\el<2\elmax$, $\forall \m$. Consequently, the sampling distribution may be written $$\label{eqn:dh_sampling_distn_2}
s(\sas) = 
1 +
\sum_{\el=2\elmax}^{\infty} \:
\summ \:
\shc{s}{\el}{\m} \: \shfarg{\el}{\m}{\sas}
\spcend .$$ The harmonic coefficients of the product of the original band-limited function and the sampling distribution $\f^s = \f \cdot s $ are given by $$\shc{\f}{\el}{\m}^s =
\sum_{\saai=0}^{2\elmax-1} \:
\sum_{\sabi=0}^{2\elmax-1} \:
\qweightdh(\saaiang) \: \f(\saisang) \: \shfargc{\el}{\m}{\saisang}
\spcend ,$$ which follows from [Eqn. (\[eqn:dh\_sampling\_distn\])]{}. Notice that these harmonic coefficients are given by the quadrature rule specified in [Eqn. (\[eqn:dh\_quad\])]{} and it simply remains to prove that the harmonic coefficients of $\f^s$ agree with those of  for the harmonic range of interest ( for $0\leq\el<\elmax$). Noting [Eqn. (\[eqn:dh\_sampling\_distn\_2\])]{}, we may write $\f^s(\sas) = 
\f(\sas) + \alpha(\sas)$, where $$\alpha(\sas) = 
\sum_{\el=2\elmax}^{\infty} \:
\summ \:
\shc{s}{\el}{\m} \:
\shfarg{\el}{\m}{\sas} \:
\sum_{\el\p=0}^{\elmax-1} \:
\sum_{\m\p=-\el\p}^{\el\p} \:
\shc{\f}{\el\p}{\m\p} \:
\shfarg{\el\p}{\m\p}{\sas} 
\spcend .$$ Since the product of two spherical harmonic functions $\shfarg{\el}{\m}{\sas} \: \shfarg{\el\p}{\m\p}{\sas}$ can be written as a sum of spherical harmonics with minimum degree $| \el - \el\p| $, [@driscoll:1994] the aliasing error $\alpha(\sas)$ contains non-zero harmonic content for $\el>\elmax$ only. Aliasing is therefore outside of the harmonic range of interest and $\shc{\f}{\el}{\m}^s = \shc{\f}{\el}{\m}$ for $0\leq\el<\elmax$, $|m|<\el$, thus proving the exact quadrature rule given by [Eqn. (\[eqn:dh\_quad\])]{}.

McEwen & Wiaux Sampling Theorem
-------------------------------

The MW sampling theorem [@mcewen:fssht] follows by a factoring of rotations [@risbo:1996] and a periodic extension in colatitude , so that the orthogonality of the complex exponentials over $[0,
2\pi)$ may be exploited. This approach encodes an implicit quadrature rule on the sphere, which can then be made explicit.

The spherical harmonics are related to the Wigner functions through [@goldberg:1967] $$\label{eqn:ssh_wigner}
\sshfarg{\el}{\m}{\sas}{{{}}} =
\sqrt{\frac{2\el+1}{4\pi} } \:
\dmatbig_{\m 0}^{\el\:\cconj}(\sab,\saa  ,0)
\spcend ,$$ where the Wigner functions form the canonical orthogonal basis on the rotation group . The Wigner functions may be decomposed as [@varshalovich:1989] $$\label{eqn:d_decomp}
\dmatbig_{\m\n}^{\el}(\euls)
= {\rm e}^{-\img \m\eula} \:
\dmatsmall_{\m\n}^\el(\eulb) \:
{\rm e}^{-\img \n\eulc}
\spcend ,$$ where the rotation group is parameterised by the Euler angles $(\euls)$, with $\eula \in [0,2\pi)$, $\eulb \in [0,\pi]$ and $\eulc
\in [0,2\pi)$. The Fourier series decomposition of the -functions is given by [@nikiforov:1991] $$\label{eqn:wigner_sum_reln}
  \dlmnb = \img^{\n-\m} \sum_{\m\p=-\el}^\el
  \dlmnhalfpi{\el}{\m\p}{\m} \:
  \dlmnhalfpi{\el}{\m\p}{\n} \:
  \exp{\img \m\p \eulb}
  \spcend ,$$ with . This expression follows from a factoring of rotations. [@risbo:1996] The Fourier series representation of  given by [Eqn. (\[eqn:wigner\_sum\_reln\])]{} allows one to write the spherical harmonic expansion of $\f(\sas)$ in terms of a Fourier series expansion of the function extended appropriately to the two-torus . Noting [Eqn. (\[eqn:ssh\_wigner\])]{} – [Eqn. (\[eqn:wigner\_sum\_reln\])]{}, the forward spherical harmonic transform may be written $$\label{eqn:flm_wig}
  \fslm = 
      \img^{\m} \nl
  \summptrunc
  \dlmnhalfpim \:
    \dlmnhalfpi{\el}{\m\p}{0} \:
  \Gsmm
  \spcend ,$$ where $$\label{eqn:Gmm}
  \Gsmm = \intsaa \: \Gsmt \: \exp{-\img \m\p \saa}$$ and $$\label{eqn:Gmt}
  \Gsmt = \intsab \: \fs(\sas) \: \exp{-\img \m \sab}
  \spcend .$$ Since [Eqn. (\[eqn:Gmt\])]{} is simply a Fourier transform, the discrete and continuous orthogonality of the complex exponentials may be exploited to evaluate this integral exactly by $$\Gsmti = \frac{2 \pi}{2\elmax-1} 
 \sumsabi \: 
 \fs(\saisang) \: 
 \exp{-\img \m \sabiang}
 \spcend ,$$ where $\sab_\sabi = 2 \pi \sabi/(2\elmax-1)$, for $\sabi =
0,\dotsc,2\elmax-2$, and $\saa_\saai =
\pi(2\saai+1)/(2\elmax-1)$, for $\saai = 0,\dotsc,\elmax-1$, giving $\Nmw = (\elmax - 1) (2 \elmax - 1) + 1 \sim 2 \elmax^2$ samples on the sphere. It remains to develop a quadrature rule to evaluate [Eqn. (\[eqn:Gmm\])]{}. This is achieved by extending  to the domain $\saa \in [0,2\pi)$ through the construction $$\begin{aligned}
 \rGsmti =  
 &\quad
 \begin{cases}
    \Gsmti \: , & \saai \in \{ 0, 1, \dotsc, \elmax-1 \} \\
   (-1)^{\m} \: 
   \Gsm(\saa_{2\elmax-2-\saai}) \: , & \saai \in \{ \elmax, \dotsc, 2\elmax-2 \}
 \end{cases}
 \spcend ,\end{aligned}$$ so that  may be expressed by a Fourier series. The factor $(-1)^\m$ is required to ensure that the symmetry in the domain $[0,2\pi)$ dictated by the inverse transform is preserved. Substituting the Fourier expansion of  into [Eqn. (\[eqn:Gmm\])]{} yields $$\Gsmm = 2 \pi \sum_{\m{\p}{\p}=-(\elmax-1)}^{\elmax-1} \: \Fsmmp \: \weight(\m{\p}{\p} - \m\p)
 \spcend ,
 \label{eqn:Gmm_convolution}$$ where the weights are given by $$\weight(\m\p)
 = \intsaa \: \exp{\img \m\p \saa}
 = 
 \begin{cases}
   \: \pm \img \pi/2, & \m\p=\pm 1\\
   \: 0, & \m\p \text{ odd},\:\m\p\neq\pm1\\
   \: 2/(1-{\m\p}^2), & \m\p \text{ even}
 \end{cases}
 \spcend ,$$ with $$\Fsmm =
 \frac{1}{2 \pi(2\elmax-1)} 
 \sumsaai \: \rGsmti \:
 \exp{- \img \m\p \saaiang}
 \spcend .$$ Since the spherical harmonic coefficients  are recovered exactly, all of the information content of the function $\f(\sas)$ is captured in the finite set of samples.

The derivation above effectively gives an implicit quadrature rule for the exact integration of a band-limited function on the sphere. This quadrature rule can be written explicitly as [@mcewen:fssht] $$\label{eqn:quadrature}
\int_\sphere \dmu{\sas} \: \fs(\sas) = 
\sum_{\saai=0}^{\elmax-1} \:
\sum_{\sabi=0}^{2\elmax-2} \:
\qweightmw(\saaiang) \:
\fs(\saisang)
\spcend ,$$ where the quadrature weights are defined by $$\qweightmw(\saaiang) = 
\frac{2\pi}{2\elmax-1} 
\Bigl[ 
\weighttrans(\saaiang) 
+ (1 - \kron{t,}{\elmax-1}) \: 
\weighttrans(\saa_{2\elmax-2-\saai})
\Bigr]
\spcend ,$$ and where $\weighttrans(\saaiang)$ is the inverse discrete Fourier transform of the reflected weights $\weight(-\m\p)$: $$\weighttrans(\saaiang) =
\frac{1}{2\elmax-1} \: \sum_{\m\p=-(\elmax-1)}^{\elmax-1} \:
\weight(-\m\p) \: \exp{\img \m\p \saaiang}
\spcend .$$

COMPRESSIVE SENSING ON THE SPHERE {#sec:compressive_sensing}
=================================

Compressive sensing on the sphere has been studied recently for signals sparse in harmonic space or in a redundant set of overcomplete dictionaries. [@abrial:2007; @rauhut:2011] However, many natural signals are sparse in measures defined in the spatial domain, such as in the magnitude of their gradient. A more efficient sampling of a band-limited signal on the sphere, as afforded by the MW sampling theorem, improves both the dimensionality and sparsity of the signal in the spatial domain, which has been shown to improve the quality of compressive sampling reconstruction.[@mcewen:css2] We review this very recent work, discussing the impact of efficient sampling on the sphere in the context of a total variation (TV) inpainting problem, after first defining the discrete TV norm on the sphere.

TV Norm on the Sphere
---------------------

The continuous TV norm on the sphere is defined by $$\| \f \|_{\rm TV} 
=
\int_\sphere \dmun \:
| \nabla \f |
\spcend ,$$ where the magnitude of the gradient of the signal $\f$ is given by $$| \: \nabla \f  \:| =
\sqrt{
\Biggl ( \frac{\partial \f}{\partial \saa} \Biggr)^2
+
\frac{1}{\sin^2\saa}\Biggl (  \frac{\partial \f}{\partial \sab} \Biggr )^2
}
\spcend .$$ In practice, however, one must consider the TV norm of the sampled signal, where the samples of $\f(\sas)$ are denoted by the concatenated vector $\vect{x}\in\complex^{\N}$, where  is the number of samples on the sphere of the chosen sampling theorem (hereafter harmonic coefficients  are also represented by a concatenated vector, denoted $\vect{\hat{x}}\in\complex^{\elmax^2}$). A discrete TV norm on the sphere is defined by approximating the continuous norm in the context of either the DH or MW sampling theorem. The integral of the continuous TV norm can be approximated using the quadrature rule corresponding to the sampling theorem on the sphere adopted: $$\label{eqn:tvnorm_cts_approx}
\| \f \|_{\rm TV} 
\simeq
\sum_{\saai=0}^{\N_{\saa}-1} \:
\sum_{\sabi=0}^{\N_{\sab}-1} \:
| \nabla \f | \:
\qweight(\saa_\saai)
\spcend ,$$ where the number of samples in $(\sas)$, given by $\N_\saa$ and $\N_\sab$ respectively, and the quadrature weights $\qweight(\saa_\saai)$ depend on the choice of sampling theorem. If $| \nabla \f |$ were band-limited, then [Eqn. (\[eqn:tvnorm\_cts\_approx\])]{} would be exact. Although this is not likely to be the case, [Eqn. (\[eqn:tvnorm\_cts\_approx\])]{} is nevertheless a reasonable approximation of the continuous TV norm. The magnitude of the gradient $|\nabla \f|$ can be approximated from the samples $\vect{x}$ using finite differences, to give a discrete TV norm on the sphere that approximates the continuous norm closely: $\|
\vect{x} \|_{\rm TV} \simeq \| \f \|_{\rm TV}$. Notice that the inclusion of the quadrature weights $\qweight(\saa_\saai)$ regularises the $\sin{\saa}$ term that arises from the definition of the gradient on the sphere, eliminating numerical instabilities that would otherwise occur.

TV Inpainting on the Sphere
---------------------------

We illustrate the impact of the number of samples of the DH and MW sampling theorems on compressive sensing reconstruction with an inpainting problem, where measurements are made in the spatial domain. A test signal sparse in its gradient is constructed from a binary Earth map, smoothed to give a signal band-limited at $\elmax=32$ (see ).[^4] The real inpainting problem $\vect{y} = 
\sensmat \vect{x} + \vect{n}$ is considered, where $\nmeas$ noisy measurements $\vect{y}\in\reals^{\nmeas}$ of the signal on the sphere $\vect{x}\in\reals^{\N}$ are made. The measurement operator $\sensmat\in\reals^{\nmeas \times \N}$ represents a random masking of the signal. The noise $\vect{n}\in\reals^\nmeas$ is assumed to be independent and identically distributed Gaussian noise, with zero mean.

The TV inpainting problem is first solved directly on the sphere: $$\label{eqn:recon_spatial}
\vect{x}^\star =
\underset{\vect{x}}{\arg \min} \:
\| \vect{x} \|_{\rm TV} \:\: \mbox{such that} \:\:
\| \vect{y} - \sensmat \vect{x}\|_2 \leq \epsilon
\spcend ,$$ where the bound $\epsilon$ is related to a residual noise level estimator. By adopting the MW sampling theorem in place of the DH sampling theorem, the dimensionality and sparsity of the signal in the spatial domain is optimised. However, no sampling theorem on the sphere reaches the optimal number of samples in the spatial domain suggested by the $\elmax^2$ dimensionality of the signal in the harmonic domain. Consequently, the dimensionality of the problem is reduced by recovering the harmonic coefficients $\vect{\hat{x}}$ directly: $$\label{eqn:recon_harmonic}
\vect{\hat{x}}^\star =
\underset{\vect{\hat{x}}}{\arg \min} \:
\| {\ensuremath{\Lambda}}\vect{\hat{x}} \|_{\rm TV} \:\: \mbox{such that} \:\:
\| \vect{y} - \sensmat {\ensuremath{\Lambda}}\vect{\hat{x}}\|_2 \leq \epsilon
\spcend ,$$ where ${\ensuremath{\Lambda}}\in\complex^{\N\times\elmax^2}$ represents the inverse spherical harmonic transform; the signal on the sphere is recovered by $\vect{x}^\star = {\ensuremath{\Lambda}}\vect{\hat{x}}^\star$. For this problem the dimensionality of the signal directly recovered $\vect{\hat{x}}$ is identical for both sampling theorems, however sparsity in the spatial domain remains superior ( fewer non-zero values) for the MW sampling theorem.

Reconstruction performance is plotted in when solving the inpainting problem in the spatial and harmonic domains, through [Eqn. (\[eqn:recon\_spatial\])]{} and [Eqn. (\[eqn:recon\_harmonic\])]{} respectively, for both sampling theorems (averaged over ten simulations of random measurement operators and independent and identically distributed Gaussian noise). Strictly speaking, compressed sensing corresponds to the measurement ratio $\nmeas/\elmax^2<1$ when considering the harmonic representation of the signal. Nevertheless, experiments are extended to $\nmeas/\elmax^2
\sim 2$, corresponding to the equivalent of complete sampling on the MW grid. When solving the inpainting problem in the spatial domain we see a large improvement in reconstruction quality for the MW sampling theorem when compared to the DH sampling theorem. This is due to the enhancement in both dimensionality and sparsity afforded by the MW sampling theorem in this setting. When solving the inpainting problem in the harmonic domain, we see a considerable improvement in reconstruction quality for each sampling theorem, since the dimensionality of the recovered signal is optimal in harmonic space. For harmonic reconstructions, the MW sampling theorem remains superior to the DH sampling theorem due to the enhancement in sparsity (but not dimensionality) that it affords in this setting. In all cases, the superior performance of the MW sampling theorem is evident. In example reconstructions are shown, where the superior quality of the MW reconstruction is again clear.

![Reconstruction performance for the DH and MW sampling theorems[]{data-label="fig:snr_vs_m"}](figures/SNRvsM){width="110mm"}

CONCLUSIONS {#sec:conclusions}
===========

Although compressive sensing states that sparse or compressible signals may be acquired with fewer samples than standard sampling theorems would suggest, the sampling theorem adopted nevertheless has an important influence on the performance of compressive sensing reconstruction. In Euclidean space, Shannon’s sampling theorem provides an optimal sampling for regular grids, leading to a unique sampling theorem. On the sphere, however, no sampling theorem is optimal, with different sampling theorems requiring a differing number of samples. The MW sampling theorem[@mcewen:fssht] has been developed only recently and achieves a more efficient sampling of the sphere than alternatives, requiring fewer than half as many samples as the canonical DH sampling theorem[@driscoll:1994], while still capturing all of the information content of a band-limited signal. A reduction by a factor of two in the number of samples between the DH and MW sampling theorems has important implications for compressive sensing on the sphere, both in terms of the dimensionality and sparsity of signals. The more efficient sampling of the MW sampling theorem has been shown to enhance the performance of compressed sensing reconstruction on the sphere, as illustrated with an inpainting problem.[@mcewen:css2]

JDM is supported by the Swiss National Science Foundation (SNSF) under grant 200021-130359. YW is supported by the Center for Biomedical Imaging (CIBM) of the Geneva and Lausanne Universities, EPFL, and the Leenaards and Louis-Jeantet foundations, and by the SNSF under grant PP00P2-123438.

[^1]: Fast algorithms have been developed to compute the forward and inverse transforms rapidly for both the DH [@driscoll:1994; @healy:2003] and MW [@mcewen:fssht] sampling theorems; these algorithms are essential to facilitate the application of these sampling theorems at high band-limits.

[^2]: Gauss-Legendre (GL) quadrature can also be used to construct an efficient sampling theorem on the sphere, with samples [@mcewen:fssht]. The MW sampling theorem nevertheless requires fewer samples and so remains more efficient, especially at low band-limits. Furthermore, it is not so straightforward to define norms describing spatial priors on the GL grid since it is not equiangular. Finally, algorithms implementing the GL sampling theorem have been shown to be limited to lower band-limits and less accurate than the algorithms implementing the MW sampling theorem [@mcewen:fssht]. Consequently, we do not focus on the GL sampling theorem any further in this article.

[^3]: The original DH sampling theorem has been revisited[@healy:2003], yielding an alternative formulation with only very minor differences and that also requires $\sim4\elmax^2$ samples.

[^4]: The original Earth topography data are taken from the Earth Gravitational Model (EGM2008) publicly released by the U.S. National Geospatial-Intelligence Agency (NGA) EGM Development Team. These data were downloaded and extracted using the tools available from Frederik Simons’ webpage: <http://www.princeton.edu/geosciences/people/simons/>.
