#### Abstract

The stereographic projection determines a bijection between the two-sphere, minus the North Pole, and the tangent plane at the South Pole. This correspondence induces a unitary map between the corresponding spaces. This map in turn leads to equivalence between the continuous wavelet transform formalisms on the plane and on the sphere. More precisely, any plane wavelet may be lifted, by inverse stereographic projection, to a wavelet on the sphere. In this work we apply this procedure to orthogonal compactly supported wavelet bases in the plane, and we get continuous, locally supported orthogonal wavelet bases on the sphere. As applications, we give three examples. In the first two examples, we perform a singularity detection, including one where other existing constructions of spherical wavelet bases fail. In the third example, we show the importance of the local support, by comparing our construction with the one based on kernels of spherical harmonics.

#### 1. Introduction

Two-dimensional wavelets are by now a standard tool in image processing, under the two concurrent approaches, the Discrete Wavelet Transform (DWT), based on the concept of multiresolution analysis, and the Continuous Wavelet Transform (CWT). While the former usually leads to wavelet bases or frames, the CWT has to be discretized for numerical implementation and produces in general only frames [1, 2].

Nowadays, many situations yield data on *spherical* surfaces or closed *sphere-like surfaces*, that is, surfaces obtained from a sphere by a smooth deformation, for instance, in Earth and Space sciences (geography, geodesy, meteorology, astronomy, cosmology, etc), in crystallography (texture analysis of crystals), in medicine (some organs are regarded as sphere-like surfaces), or in computer graphics (modelling of closed surfaces as the graph of a function defined on the sphere). So one needs a suitable analysis tool for such data. As in the flat case, the Fourier transform is a standard tool, but it amounts to an expansion in spherical harmonics, whose support is the whole sphere. Fourier analysis on the sphere is thus global and cumbersome. Therefore many different methods have been proposed to replace Fourier analysis with some sort of wavelet analysis.

In addition, some data may live on more complicated manifolds, such as a *two-sheeted hyperboloid*, in cosmology, for instance (an open expanding model of the universe), or a *paraboloid*. In optics also, data on such manifolds are essential for the treatment of omnidirectional images via the catadioptric procedure, for instance, in robotic vision. This last topic is particularly relevant for engineering purposes, because of the many applications in navigation, surveillance, and visualization. In the catadioptric image processing, a sensor overlooks a mirror, whose shape may be spherical, hyperbolic, or parabolic. However, instead of projecting the data from that mirror onto a plane, an interesting alternative consists in processing them directly on the mirror, and thus wavelets on such manifolds are needed [3]. Among the three shapes, the parabolic one is the most common (think of the headlights of a car). Now this case brings us back to the topic of this paper. Indeed it has been shown in [4] that the reconstruction of the orthographic (i.e., vertical) projection from a parabolic mirror can be computed as the inverse stereographic projection from the image plane onto the unit sphere—which is precisely the tool we are going to use in the sequel for designing wavelets on the sphere.

For an efficient wavelet analysis of signals or images, including on the sphere, the following properties are desirable (a thorough discussion of this topic may be found in [5]).

(i)*Basis*. The redundancy of the frames leads to nonunique expansions. Moreover, the existing constructions of spherical frames (see e.g., [6, 7]) are computationally heavy and often applicable only to band-limited functions. In particular, for the compression of large data sets, it is essential to have wavelet bases, not frames.(ii)

*Orthogonality*. This is the most economical method since it leads to orthogonal reconstruction matrices. The inversion of such matrices, needed, for example, in data compression, is trivial (the inverse equals the adjoint). Thus, orthogonal bases are ideal for compression, but this is not always sufficient: sparsity of reconstruction matrices is still needed in the case of large data sets. We should also mention that orthogonality of a wavelet basis can be easily achieved by a Gram-Schmidt procedure, but the locality of the support is usually lost.(iii)

*Local support*. A wavelet has

*local support*if it vanishes identically outside a small region. It is

*localized*if it is negligible outside a small region, so that it may have (small, but nonzero) “tails" there. Since these tails may spread in the process of approximation of data and spoil their good localization properties, local support is definitely preferred (see the example in Figure 4). More important, local support is crucial when working with large data sets, since it yields sparse matrices.(iv)

*Continuity, smoothness*. These properties are always desirable in approximation, but not easily achieved.

There are many constructions in literature that fulfil at least one of the properties above, so we have to restrict ourselves by mentioning only some of them. In [8–11], smooth wavelet bases were obtained by different approaches. Besides smoothness, the wavelets constructed in [11] have local support. Orthogonal, smooth and localized (but not locally supported) wavelets were constructed in [12]. In [13–15], by making use of a radial projection, one has obtained either continuous, semiorthogonal, locally supported wavelet bases, or piecewise constant orthogonal locally supported wavelet bases.

However, none of these methods so far has led to wavelet bases on the sphere which are *simultaneously* continuous (or smoother), orthogonal and locally supported. The aim of the present paper is precisely to fill this gap. The method we propose consists in lifting wavelets from the tangent plane to the sphere by inverse stereographic projection. It yields simultaneously smoothness, orthogonality, local support, vanishing moments. The disadvantage is that it could entail distortions around the North Pole . Of course, mathematically the construction presented here applies to the pointed sphere (and therefore we have to avoid the data around the North Pole ), but in practice this situation is harmless. Indeed, we can choose in a region where there are no relevant data. To give an example, European climatologists routinely put the North Pole of their spherical grid in the middle of the Pacific Ocean. Actually, most of the spherical data are given in polar coordinates not containing information at the poles.

To summarize our work, we aim at detecting and quantizing *local* singularities in data, at *small* scales. This is, from the practical point of view, the main purpose of wavelet analysis. We will come back to this point in Section 5.

#### 2. Preliminaries

We consider the unit sphere and the pointed sphere , with the parametrization

Let denote the stereographic projection from the North Pole and let , be defined as

respectively. The relations between and are and, conversely,

Then it is easy to prove the following relations between , the area element of , and , the area element of :

For simplicity, we write , since the set is of measure zero) and . Then we consider the map induced by the stereographic projection, namely,

and its inverse

It is well known that is a unitary map, hence

or, equivalently,

for all The last equality allows us to construct orthogonal bases on starting from orthogonal bases in . More precisely, we will use the fact that, if the functions are orthogonal, then the functions and will be orthogonal in .

#### 3. Multiresolution Analysis (MRA) and Wavelet Bases of

In order to fix our notations, we will briefly review in this section the standard construction of 2D orthonormal wavelet bases in the flat case, starting from a multiresolution analysis (MRA) [1].

Let be a regular matrix with the properties

(a), which is equivalent to the fact that has integer entries, (b), that is, all eigenvalues of have modulus greater than 1.A multiresolution analysis of associated to is an increasing sequence of closed subspaces with and , and satisfying the following conditions:

(1)(2)there exists a function such that the set is an orthonormal basis of .As a consequence, is an orthonormal basis for

For each , let us define the space as the orthogonal complement of into , that is, . The two-dimensional wavelets are those functions which span . One can prove (see [16]) that there exist wavelets that generate an orthonormal basis of . Therefore, is an orthonormal basis of for each , and is an orthonormal basis of

A particular case is that of tensor product wavelets, corresponding to the dilation matrix and a 1D MRA with scaling function and mother wavelet . In this case, and one gets the 2D scaling function and the three wavelets

If the one-dimensional functions and have compact support, then obviously so have and . This is the case of the well-known Daubechies wavelets *db*n that we will use in the sequel.

#### 4. Multiresolution Analysis and Orthonormal Wavelet Bases of

The construction of multiresolution analysis and wavelet bases in is based on the equality (2.10). To every function , one can associate the function as

In particular, if the functions are orthogonal, so are

Then, taking and , we obtain the spherical functions

For we define the spaces as

Using (2.10) and the unitarity of the map , it is immediate that is a closed subspace of , thus a Hilbert space. Moreover, these spaces have the following properties:

(1), (2) and , (3)the set is an orthonormal basis for We will say that a sequence of subspaces of with the properties above constitutes a *multiresolution analysis* of . Unlike the construction in [13], the functions considered here are scaled versions of the same function Moreover, in condition (3), we have here a genuine orthogonality of the wavelet basis, instead of the Riesz basis obtained in [13].

*Remark 4.1. *It should be noticed that, in the spherical MRA defined here, the scale parameter runs over the whole of , whereas in the usual constructions on , runs only from 0 to . This is due to the stereographic projection, which removes one pole, so that effectively is noncompact and thus arbitrary dilations are now allowed (see, however, Section 5).

Once the multiresolution analysis is determined, we construct the wavelet spaces in the usual manner. Let denote the orthogonal complement of the coarse space in the fine space so that

One can easily prove that, for each , is an orthogonal basis for and therefore is an orthonormal basis for .

The conclusion of the analysis may be summarized as follows.

(i)If has compact support in , then has*local support*on (indeed diam as ).(ii)An

*orthonormal*2D wavelet basis leads to an

*orthonormal*spherical wavelet basis.(iii)

*Smooth*2D wavelets lead to

*smooth*spherical wavelets. (iv)In particular, plane tensor product Daubechies wavelets lead to locally supported and orthonormal wavelets on , and so do plane tensor product Haar wavelets.(v)The decomposition and reconstruction matrices needed in the spherical case are the same as in the plane 2D case, so that the latter can be used (with existing toolboxes).

#### 5. Further Comments

The relation linking the area elements of and of implies that, near the origin of the plane, corresponding to the South Pole, the areas on the sphere and on the plane are almost the same, but near the North Pole, the ratio of the areas approach infinity. In the standard cases, the elements of a planar orthogonal wavelet basis have supports of the same size all over the plane. But the corresponding spherical wavelets obtained here will have a support that almost vanishes when they are in the neighborhood of the North Pole, which corresponds to regions far away from the origin in the plane. This might clearly lead to analytical as well as numerical problems whenever the information is nearly uniform over the whole sphere, and here we see the limitation of our method. As we have precised already at the end of Section 1, we have to avoid data around the North Pole. More precisely, our method is applicable on the whole sphere, except a spherical cap . The value of may actually depend on the type of data.

We emphasize that, in practice, wavelet transforms are useful mainly for a *local* analysis. Thus we do not pretend that our method is the best for all kinds of applications. Actually, no construction of spherical wavelets is perfect, as we have mentioned in Section 1. However, we claim that our method is more efficient than most other ones when one deals with large amounts of data, in all applications where one has to make a decomposition and a reconstruction, which implies to invert a large matrix. In the case of Daubechies wavelets, this matrix becomes orthogonal and sparse. Here we see the advantage of being able to use Daubechies wavelets rather than other 2D wavelets.

We can be more precise about the latitude effect, that is, the choice of the limiting value . When we project stereographically onto the plane spherical data situated outside of the spherical cap the projected data will be situated inside the square , where In the case of tensor product wavelets, the planar grid at level has the dimension and is taken as

Here is taken such that that is, Table 1 gives the size of the planar grid for and (and therefore of the matrix which is used in the decomposition algorithm) for different values of . For instance, in our first two examples below, we will use the values and . This gives matrices of reasonable size, yet yields good results.

#### 6. Examples

In order to illustrate the merits of our construction, we present three examples.

In the first example, we take the function where

The function and its gradient are continuous, but the second partial derivative with respect to has a discontinuity on the latitudinal circle (see Figure 1(a)). Detecting properly such a discontinuity requires a wavelet with two vanishing moments at least, so that none of the existing constructions of discrete wavelets mentioned in Section 1 could detect this discontinuity.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(g)**

**(h)**

**(i)**

Thus, following [6, 17], we take the discretized spherical CWT with the spherical wavelet , defined as in (4.1), associated to the plane wavelet

The plane wavelet (6.2) has vanishing moments of order up to 3 (here too, a simpler wavelet with less vanishing moments could not detect the discontinuity). In order to compare the two methods, we present in Figure 1 the analysis of the function for three different values of (a) ; (d) ; (g) . The second column shows, in panels (b), (e), (h), the result of the analysis with the spherical wavelet associated to the Daubechies wavelet *db*3. Finally, the third column, on panels (c), (f), and (i), does the same for the analysis by the discretized CWT method with the wavelet .

Comparing the second and the third column, it is clear that the Daubechies wavelet *db*3 lifted on the sphere by (4.1) outperforms the wavelet mentioned above. The precision is much better, in the sense that the width of the detected singular curve is narrower at all three latitudes. With both methods, the detection is still efficient very close to the North Pole, but the computation with the lifted *db*3 wavelet then requires considering matrices of large dimension, as shown in Table 1.

Looking at the third column suggests to increase the resolution of the CWT method, by choosing smaller scales. In order to test this effect, we show in Figure 2 a systematic analysis of the function by the discretized CWT method with the wavelet . Panels (a), (b), (c) and (d) present the spherical CWT at smaller and smaller scales, and , respectively. The discretization grids are those considered in [6]. From the panels (a)–(c), it appears that the discontinuity along the equator is detected properly, and the precision increases as the scale decreases. However, there is a limit: when the scale is taken below , the singularity is no more detected. For instance, for (see panel (d)), the wavelet becomes too narrow and “falls in between" the discretization points, and “ripples" appear. This effect is described in detail in [17]. Here again, the spherical Daubechies wavelet performs better.

**(a)**

**(b)**

**(c)**

**(d)**

The next question to examine is the possible distortion in latitude around the North Pole, for a signal which is no longer invariant under rotation around the -axis. Thus, in the second example, we take the function

and then , the function obtained from by performing a rotation around the axis by an angle . We apply a tensor product Haar wavelet transform to , for . Let us proceed in two steps. First, the discontinuity is well detected in the tangent plane (where the wavelet transform is computed), as shown in Figure 3(a). It is indeed a circle, with uniform intensity, as it should, in virtue of the conformal character of the stereographic projection. However, when it is lifted back onto the sphere, some artifacts appear, that can be seen in Figure 3(b). Namely, the part of the discontinuity that lies closer to the North Pole is less clearly visible than the ones that are farther away from it and there are some gaps. This is caused by the unavoidable fact that the grids on the plane and on the sphere never match, the inverse projection does not conserve areas, so that some points are lost in the lifting process.

**(a)**

**(b)**

**(c)**

**(a)**

**(b)**

**(c)**

**(d)**

On the other hand, the analysis of the same function by the discretized spherical CWT, at the smallest possible scale , shows no such effect, as seen in Figure 3(c). But here too, the singularity is detected by the discrete WT with a better precision. Note that, in the CWT case, simpler planar wavelets, with less vanishing moments, can be used, but the result is the same. One can also try to analyze the function obtained from by a rotation around the axis , but we do not expect much difference with the present example.

In the third example, we consider a data set from the texture analysis of crystals; see Figure 4(b). It consists of measurements on the sphere and its main characteristic is that the values over the sphere are constant, except for some peaks. Figure 4(c) shows the approximation at level 6 using the spherical wavelet frames constructed in [7], which are localized, but not locally supported. An example of such kernel is given in Figure 4(a). Figure 4(d) shows the approximation at level 1 using the spherical wavelet associated to *db*2. One can easily see that our wavelets are more efficient in approximating the given data set and also the algorithms are much faster than the ones for spherical harmonics.

#### 7. Conclusion

We have presented a new method for constructing spherical wavelets. Its advantage over most existing methods is that it yields locally supported orthonormal bases, which is crucial when one is confronted to very large data sets. The price to pay is that one has to avoid a small region around the North Pole, whose size may depend on the data. Within that region, serious distortions may occur. However, the examples given above show that this “forbidden" region can be very small: although we have used or , much higher values can be chosen without damage. Hence we believe this limitation does not prevent applying the method to practical situations. In any case, if the data are sufficiently well localized, one can always make a rotation that brings them around the South Pole, as we already mentioned in the Introduction. After all, one should remember that wavelet analysis is primarily a *local* analysis.

We may also mention that the method described in this paper could detect singularities of arbitrarily high order, by choosing appropriate planar wavelets. By contrast, the other existing spherical wavelets do not have enough vanishing moments and thus are not able to detect such discontinuities.

As a final remark, we may point out that the method should work for any manifold with a bijective orthogonal projection onto a fixed plane, since the latter induces a unitary map between the respective spaces. Such are, for instance, the upper sheet of the two-sheeted hyperboloid with projection onto the plane or the paraboloid with the same projection (in addition, we have seen in Section 1 that the latter case may be treated also by the stereographic projection from the sphere). It might also be extended to a local analysis on more general manifolds. A systematic discussion may be found in our recent work [18].

#### Acknowledgment

The authors thank one of the referees, whose constructive criticisms have much improved the readability of the paper.