Irène Waldspurger

Phase retrieval for wavelet transforms

Keywords: Phase retrieval, scalogram, iterative algorithms, multiscale method

This article describes a new algorithm that solves a particular phase retrieval problem, with important applications in audio processing: the reconstruction of a function from its scalogram, that is, from the modulus of its wavelet transform.

It is a multiscale iterative algorithm, that reconstructs the signal from low to high frequencies. It relies on a new reformulation of the phase retrieval problem, that involves the holomorphic extension of the wavelet transform. This reformulation allows to propagate phase information from low to high frequencies.

Numerical results, on audio and non-audio signals, show that reconstruction is precise and stable to noise. The complexity of the algorithm is linear in the size of the signal, up to logarithmic factors. It can thus be applied to large signals.

I. INTRODUCTION

The spectrogram is an ubiquitous tool in audio analysis and processing, eventually after being transformed into melfrequency cepstrum coefficients (MFCC) by an averaging along frequency bands. A very similar operator, yielding the same kind of results, is the modulus of the wavelet transform, sometimes called scalogram.

The phase of the wavelet transform (or the windowed Fourier transform in the case of the spectrogram) contains information that cannot be deduced from the single modulus, like the relative phase of two notes with different frequencies, played simultaneously. However, this information does not seem to be relevant to understand the perceptual content of audio signals [START_REF] Balan | On signal reconstruction without noisy phase[END_REF][START_REF] Risset | Exploration of timbre by analysis and synthesis[END_REF], and only the modulus is used in applications. To clearly understand which information about the audio signal is kept or lost when the phase is discarded, it is natural to consider the corresponding inverse problem: to what extent is it possible to reconstruct a function from the modulus of its wavelet transform? The study of this problem mostly begun in the early 80's [START_REF] Griffin | Signal estimation from modified short-time fourier transform[END_REF][START_REF] Nawab | Signal reconstruction from short-time fourier transform magnitude[END_REF].

On the applications side, solving this problem allows to resynthesize sounds after some transformation has been applied to their scalogram. Examples include blind source separation [START_REF] Virtanen | Monaural sound source separation by nonnegative matrix factorization with temporal continuity and sparseness criteria[END_REF] or audio texture synthesis [START_REF] Bruna | Audio texture synthesis with scattering moments[END_REF].

The reconstruction of a function from the modulus of its wavelet transform is an instance of the class of phase retrieval problems, where one aims at reconstructing an unknown signal x ∈ C n from linear measurements Ax ∈ C m , whose phase has been lost and whose modulus only is available, |Ax|. These problems are known to be difficult to solve.

Two main families of algorithms exist. The first one consists of iterative algorithms, like gradient descents or alternate projections [START_REF] Fienup | Phase retrieval algorithms: a comparison[END_REF][START_REF] Gerchberg | A practical algorithm for the determination of phase from image and diffraction plane pictures[END_REF]. In the case of the spectrogram, the oldest such algorithm is due to [START_REF] Griffin | Signal estimation from modified short-time fourier transform[END_REF]. These methods are simple to implement but, because the phase retrieval problem is non-convex, they are not guaranteed to converge towards the correct solution, and often get stuck in local minima. For measurement matrices A that are chosen at random (according to a suitable distribution), this problem can be overcome with a careful initialization of the algorithm [START_REF] Netrapalli | Phase retrieval using alternating minimization[END_REF][START_REF] Candès | Phase retrieval via wirtinger flow: theory and algorithms[END_REF][START_REF] Chen | Solving random quadratic systems of equations is nearly as easy as solving linear systems[END_REF]. However, these methods do not apply when measurements are not random. In the case of the spectrogram or scalogram, the reconstructed signals tend to present auditive artifacts. For the spectrogram, the performances can be improved by applying the algorithm window by window and not to the whole signal at the same time [START_REF] Bouvrie | An incremental algorithm for signal reconstruction from stft magnitude[END_REF]]. If additional information on the nature of the audio signal is available, it can also be taken into account in the algorithm [START_REF] Achan | Probabilistic inference of speech signals from phaseless spectrograms[END_REF][START_REF] Eldar | Sparse phase retrieval from short-time fourier measurements[END_REF]. Nevertheless, in the general case, the reconstruction results are still perfectible.

More recently, convexification methods have been proposed [START_REF] Chai | Array imaging using intensity-only measurements[END_REF][START_REF] Candès | Phaselift: exact and stable signal recovery from magnitude measurements via convex programming[END_REF]. For generic phase retrieval problems, these methods are guaranteed to return the true solution with high probability when the measurement matrix A is chosen at random. In the case of the spectrogram or scalogram, the matrix is not random and the proof does not hold. However, numerical experiments on small signals indicate that the reconstruction is in general correct [START_REF] Sun | Estimating a signal from a magnitude spectrogram via convex optimization[END_REF]Waldspurger et al., 2015]. Unfortunately, these methods have a high complexity, making them difficult to use for phase retrieval problems whose size exceeds a few hundred.

In this article, we present a new algorithm for the reconstruction of a function from its scalogram. As convexification methods, it offers a reconstruction of high quality. However, it has the complexity of an iterative method (roughly proportional to the size of the signal, up to logarithmic factors) and can be applied to large signals. The memory it requires is also proportional to the size of the signal.

The algorithm is multiscale: it performs the reconstruction frequency band by frequency band, from the lowest frequencies up to the highest ones.

The main idea of this algorithm is to introduce an equivalent formulation of the phase retrieval problem (by using the analytic extension of the wavelet transform). This reformulation has two advantages.

• It gives an explicit method to propagate towards higher scales the phase information reconstructed at lower scales. Thus, at each scale, we have a good initial guess for the subproblem that we have to solve (in the spirit of [START_REF] Netrapalli | Phase retrieval using alternating minimization[END_REF][START_REF] Candès | Phase retrieval via wirtinger flow: theory and algorithms[END_REF], but for non-random measurements).

• The local optimization algorithm naturally derived from the reformulation, although non-convex, seems very robust to the problem of local minima. In a particular case, if we allow for a small modification, we can even prove that it has no spurious local minimum.

Additionally,

• we introduce a multigrid error correction method, to detect and correct eventual errors of the reconstruction algorithm afterwards; • we use our algorithm to numerically study the intrinsic stability of the phase retrieval problem, and highlight the role played by the sparsity or non-sparsity of the wavelet coefficients.

Section II is devoted to definitions and notations. In Section III, we discuss the well-posedness of the phase retrieval problem. We explain our reformulation in Section IV, prove its equivalence with the classical formulation, and describe its advantages. In Section V, we describe the resulting algorithm. In Section VI, we discuss the superiority of multiscale algorithms over non-multiscale ones. Finally, in Section VII, we give numerical results, and empirically study the stability of the underlying phase retrieval problem.

The source code, and several reconstruction examples, are available at: http://www-math.mit.edu/ waldspur/wavelets phase retrieval.html II. DEFINITIONS AND ASSUMPTIONS All signals f [n] are of finite length N . Their discrete Fourier transform is defined by

f [k] = N -1 n=0 f [n]e -2πi kn N , k = 0, . . . , N -1.
The circular convolution of two signals f and g is

f g[k] = N -1 n=0 f [n]g[k -n], k = 0, . . . , N -1, where we set g[k -n] = g[k -n + N ] when k -n < 0.
We define a family of wavelets (ψ j ) 0≤j≤J by

ψj [k] = ψ(a j k), k = 0, ..., N -1.
where the dilation factor a can be any number in ]1; +∞[ and ψ : R → C is a fixed mother wavelet. We assume that J is sufficiently large so that ψJ is negligible outside a small set of points. An example is shown in Figure 1.

The wavelet transform is defined by:

∀f ∈ C N , W f = {f ψ j } 0≤j≤J .
The problem we consider here consists in reconstructing functions from the modulus of their wavelet transform:

Reconstruct f from {|f ψ j |} 0≤j≤J .

(1) Multiplying a function by a unitary complex does not change the modulus of its wavelet transform, so we only aim at reconstructing functions up to multiplication by a unitary complex, that is up to a global phase.

All signals are assumed to be analytic:

f [k] = 0 when N/2 < k ≤ N -1. (2)
Equivalently, we could assume the signals to be real but set the ψj [k] to zero for N/2 < k ≤ N -1.

III. WELL-POSEDNESS OF THE INVERSE PROBLEM

Developing an algorithm to solve Problem (1) makes sense only if the problem is well-posed: f must be uniquely determined from {|f ψ j |} 0≤j≤J , up to a global phase, and reconstruction must satisfy a form of stability to noise.

Unfortunately, determining whether a given phase retrieval problem is well-posed or not is a notoriously difficult question that, despite decades of works, has been solved only in a handful of cases. In the case of the scalogram (or spectrogram), a theoretical study of the well-posedness, for relatively general wavelets, seems out of reach. However, we can reasonably conjecture that the inverse problem is well-posed for generic wavelet families.

Indeed, it has been proven in several cases that reconstruction is unique, as soon as the sampling rate is high enough (for Cauchy wavelets [START_REF] Mallat | Phase retrieval for the Cauchy wavelet transform[END_REF] and for real band-limited wavelets [Alaifari et al., 2016a]; in the case of the spectrogram, it is known, for general windows, that at least almost all signals are uniquely determined [START_REF] Jaganathan | STFT phase retrieval: Uniqueness guarantees and recovery algorithms[END_REF]). By contrast, there is no case in which uniqueness is known not to hold, except for families which exhibit a severe lack of redundancy (Shannon wavelets, for example).

It has been established in [START_REF] Mallat | Phase retrieval for the Cauchy wavelet transform[END_REF]] that, when reconstruction is unique, it is also stable, in the sense that the reconstruction operator is continuous (although it may not be uniformly continuous). Moreover, in the two cases where this question of stability could be studied in more detail (the Cauchy wavelet transform in [START_REF] Mallat | Phase retrieval for the Cauchy wavelet transform[END_REF]] and the Gabor spectrogram in [Alaifari et al., 2016b]), it was shown that the reconstruction operator actually enjoyed a stronger property of local stability (the stability constants depend on the wavelet family; reconstruction is more stable when the wavelets are more redundant, that is, when the overlap between their Fourier supports is large).

In addition to these theoretical works, various numerical studies support our conjecture regarding the well-posedness of the phase retrieval problem. In [Nawab et al., 1983, IV-C], for example, it is observed that, when the reconstruction algorithm succeeds in finding a signal whose spectrogram matches the target, then this signal is indistinguishable from the correct solution. The experiments in [START_REF] Griffin | Signal estimation from modified short-time fourier transform[END_REF][START_REF] Virtanen | Monaural sound source separation by nonnegative matrix factorization with temporal continuity and sparseness criteria[END_REF] show that this is also true when the spectrogram has been subject to various deformations, indicating that reconstruction enjoys a form of stability. For short signals, stability is also numerically established in [START_REF] Jaganathan | STFT phase retrieval: Uniqueness guarantees and recovery algorithms[END_REF]. Our own experiments in Section VII are in line with these works: we never observe uniqueness issues, and we empirically demonstrate in Paragraph VII-C that reconstruction satisfies the previously mentioned local stability property.

IV. REFORMULATION OF THE PHASE RETRIEVAL PROBLEM

In the first part of this section, we reformulate the phase retrieval problem for the wavelet transform, by introducing two auxiliary wavelet families.

We then describe the two main advantages of this reformulation. First, it allows to propagate the phase information from the low-frequencies to the high ones, and so enables us to perform the reconstruction scale by scale. Second, from this reformulation, we can define a natural objective function to locally optimize approximate solutions. Although non-convex, this function has few local minima; hence, the local optimization algorithm is efficient.

A. Introduction of auxiliary wavelets and reformulation

Let us fix r ∈]0; 1[, and set r j = r a j . We define:

∀k = 0, ..., N -1, ψlow j [k] = ψj [k]r k j ; ψhigh j [k] = ψj [k]r -k j .
This definition is illustrated by Figure 2. The wavelet ψlow j has a lower characteristic frequency than ψj and ψhigh j a higher one. The following theorem explains how to rewrite a condition on the modulus of f ψ j as a condition on f ψ low j and f ψ high j .

Theorem IV.1. Let j ∈ {0, ..., J} and g j ∈ (R + ) N be fixed. Let Q j be the function whose Fourier transform is

Qj [k] = r k j g 2 j [k], (3)
∀k = N 2 -N + 1, ..., N 2 .
For any f ∈ C N satisfying the analyticity condition (2), the following two properties are equivalent:

1) |f ψ j | = g j ; 2) (f ψ low j )(f ψ high j ) = Q j .
Proof. The proof consists in showing that the second equality is the analytic extension of the first one, in a sense that will be precisely defined.

For any function h : {0, ..., N -1} → C, let P (h) be

∀z ∈ C, P (h)(z) = N 2 k= N 2 -N +1 ĥ[k]z k .
Up to a change of coordinates, P (|f ψ j | 2 ) and P (g 2 j ) are equal to P ((f ψ low j )(f ψ high j

)) and P (Q j ):

Lemma IV.2. For any f satisfying the analyticity condition (2), and for any z ∈ C:

P (|f ψ j | 2 )(r j z) = P ((f ψ low j )(f ψ high j ))(z)
and P (g 2 j )(r j z) = P (Q j )(z). This lemma is proved in the appendix A. It implies the result because then,

|f ψ j | = g j ⇐⇒ |f ψ j | 2 = g 2 j ⇐⇒ ∀z, P (|f ψ j | 2 )(z) = P (g 2 j )(z) ⇐⇒ ∀z, P (|f ψ j | 2 )(r j z) = P (g 2 j )(r j z) ⇐⇒ ∀z, P ((f ψ low j )(f ψ high j ))(z) = P (Q j )(z) ⇐⇒ (f ψ low j )(f ψ high j ) = Q j .
By applying simultaneously Theorem IV.1 to all indexes j, we can reformulate the phase retrieval problem |f ψ j | = g j , ∀j in terms of the f ψ low j 's and f ψ high j 's.

Corollary IV.3 (Reformulation of the phase retrieval problem). Let (g j ) 0≤j≤J be a family of signals in (R + ) N . For each j, let Q j be defined as in (3). Then the following two problems are equivalent.

Find f satisfying (2) such that:

∀j, |f ψ j | = g j .

⇐⇒

Find f satisfying (2) such that:

∀j, (f ψ low j )(f ψ high j ) = Q j . (4)

B. Phase propagation across scales

This new formulation yields a natural multiscale reconstruction algorithm, in which one reconstructs f frequency band by frequency band, starting from the low frequencies.

Indeed, once f ψ J , ..., f ψ j+1 have been reconstructed, it is possible to estimate f ψ low j by deconvolution. This deconvolution is stable to noise because, if r is sufficiently small, then the frequency band covered by ψ low j is almost included in the frequency range covered by ψ J , ..., ψ j+1 (see Figure 2). From f ψ low j , one can reconstruct f ψ high j , using (4):

f ψ high j = Q j f ψ low j . (5)
Finally, one reconstructs f ψ j from f ψ high j and f ψ low j . The classical formulation of the phase retrieval problem does not allow the conception of such a multiscale algorithm. Indeed, from f ψ J , ..., f ψ j+1 , it is not possible to directly estimate f ψ j : it would require performing a highly unstable deconvolution. The introduction of the two auxiliary wavelet families is essential.

C. Local optimization of approximate solutions

From the reformulation (4), we can define a natural objective function for the local optimization of approximate solutions to the phase retrieval problem. Despite its nonconvexity, this function empirically seems to have less local minima than more classical objective functions. It is thus easier to globally minimize.

The objective function has 2J + 3 variables: (h low j ) 0≤j≤J , (h high j ) 0≤j≤J and the analytic function f . The intuition is that f is the signal we aim at reconstructing and the h low j , h high j correspond to the f ψ low j 's and f ψ high j 's. The objective function is

obj(h low J , ..., h low 0 , h high J , ..., h high 0 , f ) = J j=0 ||h low j h high j -Q j || 2 2 + λ J j=0 ||f ψ low j -h low j || 2 2 + ||f ψ high j -h high j || 2 2 . ( 6
)
We additionally constrain the variables (h low j ) 0≤j≤J and (h high j ) 0≤j≤J to satisfy ∀j = 0, ..., J -1,

h low j ψ high j+1 = h high j+1 ψ low j . (7)
The first term of the objective ensures that the equalities (4) are satisfied, while the second term and the additional constraint (7) enforce the fact that the h low j 's and h high j 's must be the wavelet transforms of the same function f . 6) (solid line) and (8) (dashed line), as a function of -log(r). The second proportion is independent from r. Signals were realizations of Gaussian random processes, with 500 points; the wavelets were dyadic Morlet and there was no noise.

The parameter λ is a positive real number. In our implementation, we choose a small λ, so that the first term dominates over the second one.

If we had not introduced auxiliary wavelets, an objective function that we could have considered would have been classical obj(h J , ..., h 0 , f )

= J j=0 || |h j | 2 -g 2 j || 2 2 + λ J j=0 ||f ψ j -h j || 2 2 . (8)
However, this function seems to have more local minima, and, when minimized with a first-order method, does not produce reconstruction results of the same quality as Equation ( 6). This phenomenon is illustrated by Figure 3: starting from random initializations, we have used a first-order method to minimize one or the other objective function, and have recorded how many times we have hit a local minimum, instead of the global one. Provided that r is well chosen, we see that our new formulation performs better than the classical one.

D. Analysis for Cauchy wavelets

In this section, we give a partial theoretical justification of the fact that introducing auxiliary wavelets reduces the number of local minima. Indeed, in the case of Cauchy wavelets and in the absence of noise, we show that a modified version of our objective function (6) has no spurious critical point, in the limit λ → 0 + .

The idea is that we have an explicit expression for the manifold of stable critical points of the first part of the objective function. This expression is sufficiently simple so that we can easily compute the critical points of the (modified) second part, when restricted to that manifold.

Cauchy wavelets are parametrized by numbers p 1 , p 2 > 0:

ψ(x) = x p1 e -p2x ∀x ∈ R + .
We define (ψ low j ) 0≤j≤J , (ψ high j ) 0≤j≤J with the following choice for r: r = e -p2( a-1 a+1 ) .

With this choice, we see that

ψ high j+1 = a p1 ψ low j j = 0, . . . , J -1,
and we naturally replace Condition ( 7) by

∀j = 0, . . . , J -1, h high j+1 = a p1 h low j . (9)
In the limit λ → 0 + , the critical points of the objective function are included in the critical points of the first term

obj 1 (h low J , . . . , h low 0 , h high J , . . . , h high 0 ) = J j=0 ||h low j h high j -Q j || 2 2 .
We say that obj 1 is strict saddle at a critical point if its Hessian has at least one strictly negative eigenvalue. Strict saddle critical points can be escaped from by first-order methods [START_REF] Ge | Escaping from saddle points -online stochastic gradient for tensor decomposition[END_REF][START_REF] Lee | Gradient descent converges to minimizers[END_REF], so their presence is not too harmful.

The critical points that are not strict saddle have a simple structure, described by the following proposition (proven in Appendix B).

Proposition IV.4. Let f rec be the function that we aim at reconstructing. We assume that, for any j = 0, . . . , J, the functions f rec ψ low j and f rec ψ high j have no zero entry. Then, on the vector space defined by Equation (9), the critical points of obj 1 that are not strict saddle are

γ(f rec ψ low J ), 1 γ (f rec ψ low J-1 ), γ(f rec ψ low J-2 ), . . . , 1 γ (f rec ψ high J ), γ(f rec ψ high J-1 ), 1 γ (f rec ψ high J-2 ), . . . , such that γ ∈ (C -{0}) N def = C 1 .
The second member of the objective function ( 6) is difficult to directly analyze, but a modified version is approachable. We define

obj 2 (h low J , . . . , h high 0 , f ) = j≡J[2] ||f ψ low j -h low j || 2 2 + j ≡J[2] ||f ψ high j -h high j || 2 2 .
In the limit λ → 0 + , the minimization problem becomes:

min obj 2 (h low J , . . . , h high 0 , f ) s.t. (h low J , . . . , h high 0 , f ) ∈ C 1 × A, ( 10
)
where A is the set of signals satisfying the analyticity conditions (2).

The following theorem (proven in Appendix B), ensures that all critical points of this problem are the correct solution, up to affine scaling.

Theorem IV.5. We assume J ≥ 2, and keep the hypotheses of Proposition IV.4.

Let (h low J , . . . , h low 0 , h high J , . . . , h high 0 , f ) be a critical point of obj 2 on the real manifold C 1 × A. There exists α ∈ C - {0}, b ∈ C such that f = αf rec + b.
The modified problem ( 10) is however less stable to noise that the original version (6). In applications, we thus stick to the original formulation.

V. DESCRIPTION OF THE ALGORITHM

In this section, we describe our implementation of the multiscale reconstruction algorithm introduced in Section IV. We explain the general organization in Paragraph V-A. We then describe our exhaustive search method for solving phase retrieval problems of very small size (paragraph V-B), which our algorithm uses to initialize the multiscale reconstruction. In Paragraph V-C, we describe an additional multigrid correction step.

A. Organization of the algorithm

We start by reconstructing f ψ J from |f ψ J | and |f ψ J-1 |. We use an exhaustive search method, described in the next paragraph V-B, which takes advantage of the fact that ψJ and ψJ-1 have very small supports.

We then reconstruct the components of the wavelet transform scale by scale, as described in Section IV.

At each scale, we reconstruct f ψ low j by propagating the phase information coming from f ψ J , ..., f ψ j+1 (as explained in Paragraph IV-B). This estimation can be imprecise, so we refine it by local optimization, using the objective function defined in Paragraph IV-C, from which we drop all the terms with higher scales than j. The local optimization algorithm we use in the implementation is L-BFGS [START_REF]Updating quasi-newton matrices with limited storage[END_REF], a low-memory approximation of a second order method.

We then reconstruct f ψ high j by the equation ( 5).

At the end of the reconstruction, we run a few steps of the classical Gerchberg-Saxton algorithm to further refine the estimation.

The pseudo-code 1 summarizes the structure of the implementation.

B. Reconstruction by exhaustive search for small problems

In this paragraph, we explain how to reconstruct f ψ j from |f ψ j | and |f ψ j-1 | by exhaustive search, in the case where the support of ψj and ψj-1 is small. This is the method we use to initialize our multiscale algorithm. It is also useful for the multigrid error correction step described in the next paragraph V-C.

Lemma V.1. Let m ∈ R N and K ∈ N * be fixed. We consider the problem Find g ∈ C N s.t. |g| = m
and Supp(ĝ) ⊂ {1, ..., K}.

Algorithm 1 overview of the algorithm Input: {|f ψ j |} 0≤j≤J 1: Initialization: reconstruct f ψ J by exhaustive search.

2: for all j = J : (-1) : 0 do 3:

Estimate f ψ low j by phase propagation.

4:

Refine the values of f ψ low J , ..., f ψ low j , f ψ high J , ..., f ψ high j+1 by local optimization.

5:

Do an error correction step.

6:

Refine again.

7:

Compute f ψ high j by f ψ high j = Q j /f ψ low j . 8: end for 9: Compute f . 10: Refine f with Gerchberg-Saxton.

Output: f

This problem has at most 2 K-1 solutions, up to a global phase, and there exist a simple algorithm which, from m and K, returns the list of all possible solutions.

Proof. This lemma is a consequence of classical results about the phase retrieval problem for the Fourier transform. It can for example be derived from [START_REF] Hayes | The reconstruction of a multidimensional sequence from the phase or magnitude of its fourier transform[END_REF]. We give a proof in Appendix C.

We apply this lemma to m = |f ψ j | and |f ψ j-1 |, and construct the lists of all possible f ψ j 's and of all possible f ψ j-1 's. The true f ψ j and f ψ j-1 are the only pair in these two lists which satisfy the equality

(f ψ j ) ψ j-1 = (f ψ j-1 ) ψ j .
This solves the problem.

The number of elements in the lists is exponential in the size of the supports of ψj and ψj-1 , so this algorithm has a prohibitive complexity when the supports become large. Otherwise, our numerical experiments show that it works well.

C. Error correction

When the modulus are noisy, there can be errors during the phase propagation step. The local optimization generally corrects them, if run for a sufficient amount of time, but, for the case where some errors are left, we add, at each scale, a multigrid error correction step. This step does not totally remove the errors but greatly reduces their amplitude.

1) Principle: First, we determine the values of n for which f ψ low j [n] and f ψ high j+1 [n] seem to have been incorrectly reconstructed. We use the fact that f ψ low j and f ψ high j+1 must satisfy

(f ψ low j ) ψ high j+1 = (f ψ high j+1 ) ψ low j .
The points where this equality does not hold provide a good estimation of the places where the values of f ψ low j and f ψ high j+1 are erroneous.

We then construct a set of smooth "windows" w 1 , ..., w S , whose supports cover the interval on which errors have been found (see figure 4), such that each window has a small support. For each s, we reconstruct (f ψ low j ).w s and (f ψ high j+1 ).w s , by expressing these functions as the solutions to phase retrieval problems of small size, which we can solve by the exhaustive search method described in Paragraph V-B.

If w s is sufficiently smooth so that it can be considered constant on the support of ψ j , we have, for k = 0, . . . , N -1,

|(f.w s ) ψ j [k]| = N -1 n=0 ψ j [n]w s [k -n]f [k -n] ≈ N -1 n=0 ψ j [n]w s [k]f [k -n] = w s [k] |f ψ j [k]| .
The same reasoning applies for ψ j+1 :

|(f.w s ) ψ j | ≈ w s |f ψ j |; |(f.w s ) ψ j+1 | ≈ w s |f ψ j+1 |.
The wavelets ψ j and ψ j+1 have a small support in the Fourier domain, if we truncate them to the support of w s , so we can solve this problem by exhaustive search, and reconstruct (f.w s ) ψ j and (f.w s ) ψ j+1 .

From (f.w s ) ψ j and (f.w s ) ψ j+1 , we reconstruct (f ψ low j ).w s ≈ (f.w s ) ψ low j and (f ψ high j+1 ).w s ≈ (f.w s ) ψ high j+1 by deconvolution.

There is a trade-off to reach in the size of the support of the w s : the approximation error in |(f.w s ) ψ j | ≈ w s |f ψ j | is inversely proportional to this size, but the complexity of the exhaustive search becomes prohibitive when the size is too large. In our implementation, we adjust it so that the Fourier support of ψ j (truncated to the support of w s ) is of size 12.

2) Usefulness of the error correction step: The error correction step does not perfectly correct the errors, but greatly reduces the amplitude of large ones. Figure 5 shows an example of this phenomenon. It deals with the reconstruction of a difficult audio signal, representing a human voice saying "I'm sorry". Figure 5a shows f ψ low 7 after the multiscale reconstruction at scale 7, but before the error correction step. The reconstruction presents large errors. So the error correction step must be used when large errors are susceptible to occur, and turned off otherwise: it makes the algorithm faster without reducing its precision.

Figure 6 illustrates this affirmation by showing the mean reconstruction error for the same audio signal as previously. When 200 iterations only are allowed at each local optimization step, there are large errors in the multiscale reconstruction; the error correction step significantly reduces the reconstruction error. When 2000 iterations are allowed, all the large errors can be corrected during the local optimization steps and the error correction step is not useful.

VI. MULTISCALE VERSUS NON-MULTISCALE

Our reconstruction algorithm has very good reconstruction performances, mainly because it uses the reformulation of the phase retrieval problem introduced in Section IV. However, the quality of its results is also due to its multiscale structure. It is indeed known that, for the reconstruction of functions from their spectrogram or scalogram, multiscale algorithms perform better than non-multiscale ones [START_REF] Bouvrie | An incremental algorithm for signal reconstruction from stft magnitude[END_REF][START_REF] Bruna | Scattering representations for recognition[END_REF]. In this section, we propose two justifications for this phenomenon (paragraph VI-A). We then introduce a multiscale version of the classical Gerchberg-Saxton algorithm, and numerically verify that it yields better reconstruction results than the usual non-multiscale version (paragraph VI-B).

A. Advantages of the multiscale reconstruction

At least two factors can explain the superiority of multiscale methods, where the f ψ j 's are reconstructed one by one, and not all at the same time.

First, they can partially remedy the possible ill-conditioning of the problem. In particular, if the f ψ j 's have very different norms, then a non-multiscale algorithm will be more sensitive to the components with a high norm. It may neglect the information given by |f ψ j |, for the values of j such that this function has a small norm. With a multiscale algorithm where all the |f ψ j |'s are successively considered, this happens less frequently.

Second, iterative algorithms, like Gerchberg-Saxton, are very sensitive to the choice of their starting point (hence the care given to their initialization in the literature [START_REF] Netrapalli | Phase retrieval using alternating minimization[END_REF][START_REF] Candès | Phase retrieval via wirtinger flow: theory and algorithms[END_REF]). If all the components are reconstructed at the same time and the starting point is randomly chosen, the algorithm almost never converges towards the correct solution: it gets stuck in a local minima.

In a multiscale algorithm, the starting point at each scale can be chosen so as to be consistent with the values reconstructed at lower scales; it yields much better results.

B. Multiscale Gerchberg-Saxton

To justify the efficiency of the multiscale approach, we introduce a multiscale version of the classical Gerchberg-Saxton algorithm [START_REF] Gerchberg | A practical algorithm for the determination of phase from image and diffraction plane pictures[END_REF] (by alternate projections) and compare its performances with the nonmultiscale algorithm.

The multiscale algorithm reconstructs f ψ J by exhaustive search (paragraph V-B).

Then, for each j, once f ψ J , ..., f ψ j+1 are reconstructed, an initial guess for f ψ j is computed by deconvolution. The frequencies of f ψ j for which the deconvolution is too unstable are set to zero. The regular Gerchberg-Saxton algorithm is then simultaneously applied to f ψ J , ..., f ψ j .

We test this algorithm on realizations of Gaussian random processes (see Section VII-B for details), of various lengths. On Figure 7, we plot the mean reconstruction error obtained with the regular Gerchberg-Saxton algorithm and the error obtained with the multiscale version (see Paragraph VII-A for the definition of the reconstruction error).

None of the algorithms is able to perfectly reconstruct the signals, in particular when their size increases. However, the multiscale algorithm clearly yields better results, with a mean error approximately twice smaller.

VII. NUMERICAL RESULTS

In this section, we describe the behavior of our algorithm. We compare it with Gerchberg-Saxton and with PhaseLift. We show that it is much more precise than Gerchberg-Saxton. It is comparable with PhaseLift in terms of precision, but significantly faster, so it allows to reconstruct larger signals.

The performances strongly depend on the type of signals we consider. The main source of difficulty for our algorithm is the presence of small values in the wavelet transform, especially in the low frequencies.

Indeed, the reconstruction of f ψ high j by Equation ( 5) involves a division by f ψ low j . When f ψ low j has small values, this operation is unstable and induces errors.

As we will see in Section VII-C, the signals whose wavelet transform has many small values are also the signals for which the phase retrieval problem is the least stable (in the sense that two functions can have wavelet transforms almost equal in modulus without being close in l 2 -norm). This suggests that this class of functions is intrinsically the most difficult to reconstruct; it is not an artifact of our algorithm.

We describe our experimental setting in Paragraph VII-A. In Paragraph VII-B, we give detailed numerical results for various types of signals. In Paragraph VII-C, we use our algorithm to investigate the stability to noise of the underlying phase retrieval problem. Finally, in Paragraph VII-D, we study the influence of various parameters on the quality of the reconstruction.

A. Experimental setting

At each reconstruction trial, we choose a signal f and compute its wavelet transform {|f ψ j |} 0≤j≤J . We corrupt it with a random noise n j :

h j = |f ψ j | + n j . (11)
We measure the amplitude of the noise in l 2 -norm, relatively to the l 2 -norm of the wavelet transform:

amount of noise = j ||n j || 2 2 j ||f ψ j || 2 2 . ( 12
)
In all our experiments, the n j are realizations of Gaussian white noises. We run the algorithm on the noisy wavelet transform {h j } 0≤j≤J . It returns a reconstructed signal f rec . We quantify the reconstruction error by the difference, in relative l 2 -norm, between the modulus of the wavelet transform of the original signal f and the modulus of the wavelet transform of the reconstructed signal f rec :

reconstruction error = j || |f ψ j | -|f rec ψ j | || 2 2 j ||f ψ j || 2 2 . (13)
Alternatively, we could measure the difference between f and f rec (up to a global phase):

error on the signal = inf φ∈R ||e iφ f -f rec || 2 ||f || 2 . ( 14
)
But we know that the reconstruction of a function from the modulus of its wavelet transform is not stable to noise [START_REF] Mallat | Phase retrieval for the Cauchy wavelet transform[END_REF]]. So we do not hope the difference between f and f rec to be small. We just want the algorithm to reconstruct a signal f rec whose wavelet transform is close to the wavelet transform of f , in modulus. Thus, the reconstruction error ( 13) is more relevant to measure the performances of the algorithm.

In all the experiments, unless otherwise specified, we use dyadic Morlet wavelets, to which we subtract Gaussian functions of small amplitude so that they have zero mean:

ψ(ω) = exp(-p(ω -1) 2 ) -β exp(-pω 2 )
where β > 0 is chosen so that ψ(0) = 0 and the parameter p is arbitrary (it controls the frequency bandwidth of the wavelets). For N = 256, our family of wavelets contains eight elements, which are plotted on Figure 19a. The performances of the algorithm strongly depend on the choice of the wavelet family; this is discussed in Paragraph VII-D1.

The maximal number of iterations per local optimization step is set to 10000 (with an additional stopping criterion, so that the 10000-th iteration is not always reached). We study the influence of this parameter in Paragraph VII-D2.

The error correction step described in Paragraph V-C is always turned on.

Gerchberg-Saxton is applied in a multiscale fashion, as described in Paragraph VI-B, which yields better results than the regular implementation.

We use PhaseLift [START_REF] Candès | Phase retrieval via matrix completion[END_REF] with ten steps of reweighting, followed by 2000 iterations of the Gerchberg-Saxton algorithm. In our experiments with PhaseLift, we only consider signals of size N = 256. Handling larger signals is difficult with a straightforward Matlab implementation.

B. Results

We describe four classes of signals, whose wavelet transforms have more or less small values. For each class, we plot the reconstruction error of our algorithm, Gerchberg-Saxton and PhaseLift as a function of the noise error.

1) Realizations of Gaussian random processes: We first consider realizations of Gaussian random processes. A signal f in this class is defined by

f [k] = X k √ k + 1 if k ∈ {1, ..., N/2}; = 0 if not.
where X 1 , ..., X N/2 are independent realizations of complex Gaussian centered variables. The role of the √ k + 1 is to ensure that all components of the wavelet transform approximately have the same l 2 -norm (in expectation). An example is displayed on Figure 8, along with the modulus of its wavelet transform.

The wavelet transforms of these signals have few small values, disposed in a seemingly random pattern. This is the most favorable class for our algorithm.

The reconstruction results are shown in Figure 9. Even for large signals (N = 10000), the mean reconstruction error is proportional to the input noise (generally 2 or 3 times smaller); this is the best possible result. The performances of PhaseLift are exactly the same, but Gerchberg-Saxton often fails.

2) Lines from images: The second class consists in lines randomly extracted from photographs. These signals have oscillating parts (corresponding to the texture zones of the initial image) and smooth parts, with large discontinuities in between. Their wavelet transforms generally contain a lot a small values, but, as can be seen in Figure 10, the distribution of these small values is particular. They are more numerous at high frequencies and the non-small values tend to concentrate on vertical lines of the time-frequency plane. This distribution is favorable to our algorithm: small values in the wavelet transform are mostly a problem when they are in the low frequencies and prevent the correct initialization of the reconstruction at medium or high frequencies. Small values at high frequencies are not a problem. Indeed, as in the case of Gaussian signals, the reconstruction error is proportional to the input noise (figure 11). This is also the case for PhaseLift but not for Gerchberg-Saxton.

3) Sums of a few sinusoids: The next class of signals contains sums of a few sinusoids, multiplied by a window function w to avoid boundary effects. Formally, a signal in this class is of the form

f [n] =   N/2 k=1 α k exp i 2πkn N   × w[n].
where the α k are zero with high probability and realizations of complex Gaussian centered variables with small probability.

The wavelet transforms of these signals often have components of very small amplitude, which may be located at any frequential scale (figure 12). This can prevent the reconstruction.

The results are on Figure 13. Our algorithm performs much better than Gerchberg-Saxton but the results are not as good as for the two previous classes of signals.

In most reconstruction trials, the signal is correctly reconstructed, up to an error proportional to the noise. But, with a small probability, the reconstruction fails. The same phenomenon occurs for PhaseLift.

The probability of failure seems a bit higher for PhaseLift than for our algorithm. For example, when the signals are of size 256 and the noise has a relative norm of 0.01%, the reconstruction error is larger than the noise error 20% of the time for PhaseLift and only 10% of the time for our algorithm. However, PhaseLift has a smaller mean reconstruction error because, in these failure cases, the result it returns, although not perfect, is more often close to the truth: the mean reconstruction error in the failure cases is 0.2% for PhaseLift versus 1.7% for our algorithm.

4) Audio signals: Finally, we test our algorithm on real audio signals. These signals are difficult to reconstruct because they do not contain very low frequencies (as the human ear cannot hear them, these frequencies are not included in the recordings), so the first components of their wavelet transforms are very small.

The reconstruction results may vary from one audio signal to the other. We focus here on two representative examples.

The first signal is an extract of five seconds of a musical piece played by an orchestra (the Flight of the Bumblebee, by Rimsky-Korsakov). Figure 14a shows the modulus of its wavelet transform. It has 16 components and 9 of them (the ones with lower characteristic frequencies) seem negligible, compared to the other ones. However, its non-negligible components have a moderate number of small values.

The second signal is a recording of a human voice saying "I'm sorry" (figure 14b). The low-frequency components of its wavelet transform are also negligible, but even the highfrequency components tend to have small values, which makes the reconstruction even more difficult.

The results are presented in Figures 15 and16. For relatively high levels of noise (0.5% or higher), the results, in the sense of the l 2 -norm, are satisfying: the reconstruction error is smaller or equal to the amount of noise. In the high precision regime (that is, for 0.1% of noise or less), the lack of low frequencies does not allow a perfect reconstruction. Nevertheless, the results are still good: the reconstruction error is of the order of 0.1% or 0.2% when the noise error is below 0.1%. More iterations in the optimization steps can further reduce this error. By comparison, the reconstruction error with Gerchberg-Saxton is always several percent, even when the noise is small.

C. Stability of the reconstruction

In this section, we use our reconstruction algorithm to investigate the stability of the reconstruction. From [START_REF] Mallat | Phase retrieval for the Cauchy wavelet transform[END_REF], we know that the reconstruction operator is not uniformly continuous: the reconstruction error (13) can be small (the modulus of the wavelet transform is almost exactly reconstructed), even if the error on the signal ( 14) is not small (the difference between the initial signal and its reconstruction is large).

We show that this phenomenon can occur for all classes of signals, but is all the more frequent when the wavelet transform has a lot of small values, especially in the low frequency components.

We also experimentally show that, when this phenomenon happens, the original and reconstructed signals have their wavelet transforms {f ψ j (t)} j∈Z,t∈R equal up to multiplication by a phase {e iφj (t) } j∈Z,t∈R , which varies slowly in both j and t, except maybe at the points where f ψ j (t) is close to zero. This has been conjectured in [START_REF] Mallat | Phase retrieval for the Cauchy wavelet transform[END_REF]. It is a form of "local stability": the signal and its reconstruction are close, up to a global phase, in a neighborhood of each point of the time-frequency plane.

We perform a large number of reconstruction trials, with various reconstruction parameters. This gives us a large number of pairs (f,

f rec ), such that ∀j, t, |f ψ j (t)| ≈ |f rec ψ j (t)|.
For each one of these pairs, we compute error on the modulus

= j || |f ψ j | -|f rec ψ j | || 2 2 j ||f ψ j || 2 2 ; (13) error on the signal = ||f -f rec || 2 ||f || 2 . ( 14
)
The results are plotted on Figure 17, where each point corresponds to one reconstruction trial. The x-coordinate represents the error on the modulus and the y-coordinate the error on the signal.

We always have: error on the modulus ≤ C × (error on the function).

with C a constant of the order of 1. This is not surprising because the modulus of the wavelet transform is a Lipschitz operator, with a constant close to 1.

As expected, the converse inequality is not true: the error on the function can be significantly larger than the error on the modulus. For each class, an important number of reconstruction trials yield errors such that: error on the signal ≈ 30 × error on the modulus.

For realizations of Gaussian random processes or for lines extracted from images (figures 17a and 17b), the ratio between the two errors never exceeds 30 (except for one outlier). But for sums of a few sinusoids (17c) or audio signals (17d), we may even have: error on the signal ≥ 100 × error on the modulus. These results have a theoretical justification. [START_REF] Mallat | Phase retrieval for the Cauchy wavelet transform[END_REF] explain how, from any signal f , it is possible to construct g such that |f ψ j | ≈ |g ψ j | for all j but f ≈ g in the l 2 -norm sense.

The principle of the construction is to multiply each f ψ j (t) by a phase e iφj (t) . The function (j, t) → e iφj (t) must be chosen so that it varies slowly in both j and t, except maybe at the points (j, t) where f ψ j (t) is small. Then there exist a signal g such that (f ψ j (t))e iφj (t) ≈ g ψ j (t) for any j, t. Taking the modulus of this approximate equality yields

∀j, t |f ψ j (t)| ≈ |g ψ j (t)|.
However, we may not have f ≈ g. This construction works for any signal f (unless the wavelet transform is very localized in the time frequency domain), but the number of possible {e iφj (t) } j,t is larger when the wavelet transform of f has a lot of small values, because the constraint of slow variation is relaxed at the points where the wavelet transform is small (especially when the small values are in the low frequencies). This is probably why instabilities occur for all kinds of signals, but more frequently when the wavelet transforms have a lot of zeroes.

From our experiments, it seems that the previous construction describes all the instabilities: when the wavelet transforms of f and f rec have almost the same modulus and f is not close to f rec , then the wavelet transforms of f and f rec are equal up to slow-varying phases {e iφj (t) } j,t . Figure 18. (a) modulus of the wavelet transform of a signal (b) phase difference between the wavelet transform of this signal and of its reconstruction (black points correspond to places where the modulus is too small for the phase to be meaningful) Figure 18 shows an example. The signal is a sum of sinusoids. The relative difference between the modulus is 0.3%, but the difference between the initial and reconstructed signals is more than a hundred times larger; it is 46%. The second subfigure shows the difference between the phases of the two wavelet transforms. It indeed varies slowly, in both time and frequency (actually, it is almost constant along the frequency axis), and a bit faster at the extremities, where the wavelet transform is closer to zero.

D. Influence of the parameters

In this paragraph, we analyze the importance of the two main parameters of the algorithm: the choice of the wavelets (paragraph VII-D1) and the number of iterations allowed per local optimization step (paragraph VII-D2).

1) Choice of the wavelets: Two properties of the wavelets are especially important: the exponential decay of the wavelets in the Fourier domain (so that the Q j 's (3) are correctly computed) and the amount of overlap between two neighboring wavelets (if the overlap is too small, then f ψ J , ..., f ψ j+1 contain not much information about f ψ j and the multiscale approach is less efficient). We compare the reconstruction results for four families of wavelets.

The first family (figure 19a) is the one we used in all the previous experiments. It contains dyadic Morlet wavelets. The second family (figure 19b) also contains Morlet wavelets, with a smaller bandwidth (Q-factor ≈ 8) and a dilation factor of 2 1/8 instead of 2. This is the kind of wavelets used in audio processing. The third family (figure 19c) consists in dyadic Laplacian wavelets ψ(ω) = ω 2 e 1-ω 2

. Finally, the wavelets of the fourth family (figure 19d) are (derivatives of) Gammatone wavelets.

Figure 20 displays the mean reconstruction error as a function of the noise, for two classes of signals: lines randomly extracted from natural images and audio signals.

Morlet wavelets have a fast decay and consecutive wavelets overlap well. This does not depend upon the dilation factor so the reconstruction performances are similar for the two Morlet families (figures 20a and 20b).

Laplacian wavelets are similar, but the overlap between consecutive wavelets is not as good. So Laplacian wavelets globally have the same behavior as Morlet wavelets but require significantly more computational effort to reach the same precision. Figures 20a and20b have been generated with a maximal number of iterations per optimization step equal to 30000 (instead of 10000) and the reconstruction error is still larger.

The decay of Gammatone wavelets is polynomial instead of exponential. The products Q j cannot be efficiently estimated and our method performs significantly worse. In the case of lines extracted from images (20a), the reconstruction error stagnates at 0.1%, even when the noise is of the order of 0.01%. For audio signals (20b), it is around 1% for any amount of noise. Figure 21 represents, for an audio signal, the reconstruction error as a function of this number of iterations. As the objective functions are not convex, there are no guarantees on the speed of the decay when the number of iterations increases. It can be slow and even non-monotonic. Nevertheless, it clearly globally decays.

The execution time is roughly proportional to the number of iterations. It is thus important to adapt this number to the desired application, so as to reach the necessary precision level without making the algorithm exaggeratedly slow.

VIII. CONCLUSION

We have presented an new iterative algorithm that reconstructs a signal from its scalogram. This algorithm is based on a new reformulation of the reconstruction problem, using the analytic extension of the wavelet transform. It is precise and stable to noise, and has a sufficiently low complexity to be applied to audio signals.

In future works, we plan to investigate further ways to speed up the reconstruction, including parallelization, and to test our algorithm on concrete applications, like source separation. Lemma (IV.2). For any f satisfying the analyticity condition (2),

∀z ∈ C P (|f ψ j | 2 )(r j z) = P ((f ψ low j )(f ψ high j ))(z)
and P (g 2 j )(r j z) = P (Q j )(z).

Proof. Recall that, by definition, for any

h ∈ C N , ∀z ∈ C P (h)(z) = N 2 k= N 2 -N +1 ĥ[k]z k .
So for any two signals h, H, the condition P (h)(r j z) = P (H)(z), ∀z ∈ C is equivalent to

∀k = N 2 -N + 1, ..., N 2 , ĥ[k]r k j = Ĥ[k]. (15)
Applied to h = g 2 j and H = Q j , this property yields the equality P (g 2 j )(r j z) = P (Q j )(z), ∀z ∈ C: by the definition of Q j in (3), Equation ( 15) is clearly satisfied.

Let us now show that

P (|f ψ j | 2 )(r j z) = P ((f ψ low j )(f ψ high j ))(z), ∀z ∈ C.
It suffices to prove that (15) holds, that is,

∀k ∈ N 2 -N + 1, ..., N 2 , |f ψ j | 2 [k]r k j = (f ψ low j )(f ψ high j )[k].
Indeed, because the analyticity condition (2) holds, we have for all k,

|f ψ j | 2 [k] = 1 N f ψ j f ψ j [k] = 1 N N/2 l=1 f [l] ψj [l] f [l -k] ψj [l -k] = r -k j N N/2 l=1 f [l] ψlow j [l] f [l -k] ψhigh j [l -k] = r -k j N f ψ low j f ψ high j [k] = r -k j (f ψ low j )(f ψ high j ) [k].
APPENDIX B PROOFS OF PROPOSITION IV.4 AND THEOREM IV.5 Proposition (IV.4). Let f rec be the function that we aim at reconstructing. We assume that, for any j = 0, . . . , J, the functions f rec ψ low j and f rec ψ high j have no zero entry. Then, on the vector space defined by Equation (9), the critical points of obj 1 that are not strict saddle are

γ(f rec ψ low J ), 1 γ (f rec ψ low J-1 ), γ(f rec ψ low J-2 ), . . . , 1 γ (f rec ψ high J ), γ(f rec ψ high J-1 ), 1 γ (f rec ψ high J-2 ), . . . , such that γ ∈ (C -{0}) N def = C 1 .
Proof. We first note that, by Theorem IV.1,

Q j = (f rec ψ low j )(f rec ψ high j ) j = 0, . . . , J. (16)
Let H = (h low J , . . . , h high 0 ) be a critical point that is not strict saddle. The derivative of obj 1 at H is

= h high j (h low j h high j -Q j ) g high j = h low j (h low j h high j -Q j ).
As H is a critical point on the vector space defined by constraints (9), we must have

g low J = 0; ∀j = 0, . . . , J -1, a p1 g high j+1 + g low j = 0; (17)
g high 0 = 0.
We set γ = h low J (f rec ψ low J ) -1 . Let us show by decreasing induction over j = J, . . . , 0 that h low j , h high j satisfy

h low j = γ(f rec ψ low j ) (18a)
and

h high j = 1 γ (f rec ψ high j ) if j ≡ J[2]; (18b)
h low j = 1 γ (f rec ψ low j ) (18c)
and

h high j = γ(f rec ψ high j ) if j ≡ J[2]. (18d)
We start with j = J. If we can show that h high J has no zero entry, then the equality g low J = 0 implies, by the definition of g low J :

h low J h high J = Q J = (f rec ψ low J )(f rec ψ high J ).
As f rec ψ low J and f rec ψ high J have no zero entry (by assumption), from the previous equation and the definition of γ,

h low J = γ(f rec ψ low J ) and h high J = 1 γ (f rec ψ high J ).
So let us just show that h high J has no zero entry. By contradiction, we assume that h high J

[k] = 0 for some k = 0, . . . , N -1. Let S 1 ∈ (C N ) 2(J+1) be defined as S 1 = (δ k , 0, . . . , 0, 0, . . . , 0),

where δ k ∈ C N is the vector such that δ k [l] = 0 if l = k and δ k [k] = 1. Similarly, we define S 2 = (0, δ k , . . . , 0, a p1 δ k , . . . , 0) (where a p1 δ k is the J + 2-th component).
Using the equality h high J

[k] = 0, we see that

d 2 obj 1 (H).(S 1 , S 1 ) = 0; (19) d 2 obj 1 (H).(iS 1 , iS 1 ) = 0.
However,

d 2 obj 1 (H).(S 1 , S 2 ) = -2a p1 Re (Q J [k]); d 2 obj 1 (H).(iS 1 , S 2 ) = -2a p1 Im (Q J [k]). As Q J [k] = (f ψ low J [k])(f ψ high J [k]) = 0, it means that d 2 obj 1 (H).(S 1 , S 2 ) = 0 or d 2 obj 1 (H).(iS 1 , S 2 ) = 0.
Combined with Equation ( 19), it implies that d 2 obj 1 (H) has a strictly negative eigenvalue, so H is strict saddle. This is a contradiction, and concludes the initialization.

We now assume the inductive hypothesis to hold from j + 1 to J, for some j < J, and prove it for j. From Equations ( 16), (18) and the definition of g j+1 , g high j+1 = 0. From Equation (17), g low j = 0, so, if h high j has no zero entry, we must have

h low j h high j = Q j .
Constraint ( 9) is satisfied, so, if j ≡ J[2], by Equation ( 18),

h low j = a -p1 h high j+1 = a -p1 γ(f rec ψ high j+1 ) = γ(f rec ψ low j ).
Combined with the previous equation and Equation (16), this relation also yields

h high j = 1 γ (f rec ψ high j ).
A similar reasoning holds if j ≡ J[2], so the induction hypothesis holds at rank j, provided that we show that h high j has no zero entry.

To prove this last affirmation, we reason by contradiction, and assume h high j

[k] = 0 for some k = 0, . . . , N -1. Multiplying everything by a suitable global phase, we can assume Re (h low j [k]Q j [k]) = 0 (this number is non zero, by the assumptions over f rec and by the inductive hypothesis applied to h low j = a -p1 h high j+1 ). We define

S 1 = (h low J [k]δ k , -h low J-1 [k]δ k , . . . , ±h low j [k]δ k , 0, . . . , 0, -h high J [k]δ k , h high J-1 [k]δ k , . . . , ±h high j+1 [k]δ k , 0, . . . , 0
); S 2 = (0, . . . , δ k , . . . , 0, 0, . . . , a p1 δ k , . . . , 0). (In the second definition, the non-zero components are at positions J -j + 2 and 2J -j + 2.)

We have

d 2 obj 1 (H).(S 1 , S 1 ) = 0; d 2 obj 1 (H).(S 1 , S 2 ) = -2(-1) J-j a p1 Re (h low j [k]Q j [k]) = 0. So d 2 obj 1 (H)
has at least a strictly negative eigenvalue; this is a contradiction.

This concludes the induction, and proves that all critical points that are not strict saddle have the desired form. Conversely, any point that has this form is a global minimum of obj 1 , by Equation ( 16), so it is a critical point, and is not strict saddle.

Theorem (IV.5). We assume J ≥ 2, and keep the hypotheses of Proposition IV.4.

Let (h low J , . . . , h low 0 , h high Proof. The tangent to C 1 × A at a point H = (h low J , . . . , h high J , . . . , f ) is the real vector space δh low J , -δh low J-1 , δh low J-2 , . . . , -δh high J , δh high J-1 , . . . , g , such that δ ∈ C N , g ∈ A .

If H is a critical point of obj 2 , we evaluate dobj 2 (H) along the vector of the tangent space that corresponds to δ = 1 and g = f , and we obtain 2obj 2 (H) = 0.

In particular, h low J = f ψ low J ; h low J-2 = f ψ low J-2 . So, from the definition of C 1 , there exists γ ∈ (C -{0}) N such that γ(f rec ψ low J ) = f ψ low J ;

(20a)

γ(f rec ψ low J-2 ) = f ψ low J-2 . ( 20b
)
From the definition of Cauchy wavelets and r, we deduce that ψlow J [k] = a Jp1 k p1 exp -ka J p 2 2a a + 1 , ∀k = 0, . . . , N -1.

So if we define

P (f )(X) = 1 N [N/2] l=1 l p1 f [l]X l ,
we see that, for any k = 0, . . . , N -1,

f ψ low J [k] = 1 N N -1 l=0 f [l] ψJ [l] e 2πik N l = a Jp1 P (f ) R J e 2πik N
, where R J = exp -a J p 2 2a a + 1 .

Similarly, f ψ low J-2 [k] = a (J-2)p1 P (f )(R J-2 e 2πik N ), where R J-2 = exp -a J-2 p 2 2a a + 1 .

In the same way, we have

f rec ψ low J [k] = a Jp1 P (f rec )(R J e 2πik
N );

f rec ψ low J-2 [k] = a (J-2)p1 P (f rec )(R J-2 e 2πik N ),
if we define

P (f rec )(X) = 1 N [N/2] l=1 l p1 frec [l]X l .
A consequence of these previous equalities and Equations ( 20a) and ( 20b) is that

P (f )(R J x) P (f rec )(R J x) = P (f )(R J-2 x) P (f rec )(R J-2 x) ∀x = 1, e 2πi
N , e 4πi N , . . . .

(The denominator of these fractions is never zero, because f rec ψ low J and f rec ψ low J-2 have no zero entry.) The rational fraction P (f )(X) P (frec)(X) = P (f )(X)/X P (frec)(X)/X is the quotient of two polynomials with degree at most [N/2] -1, so it is uniquely determined by its values in N complex points, and Equation (21) implies:

P (f )(x) P (f rec )(x) = P (f )(xR J-2 /R J ) P (f rec )(xR J-2 /R J ) , ∀x ∈ C.
Applying recursively this equality yields, for any s ∈ N,

P (f )(x) P (f rec )(x) = P (f ) (x(R J-2 /R J ) s ) P (f rec ) (x(R J-2 /R J ) s ) , ∀x ∈ C.
Depending on the exact degrees of P (f ) and P (f rec ), the right-hand side of this equalities converges (uniformly on every compact set) to either ∞ or a constant function, when s → +∞. The convergence to ∞ is actually impossible, because P (f ) P (frec) has finite values at almost any point of C, so P This problem has at most 2 K-1 solutions, up to a global phase, and there exist a simple algorithm which, from m and K, returns the list of all possible solutions.

Proof. We define a polynomial P (g) by

P (g)(X) = ĝ[1] + ĝ[2]X + ... + ĝ[K]X K-1 .
We show that the constraint |g| = m amounts to knowing P (g)(X)P (g)(1/X). This is in turn equivalent to knowing the roots of P (g) (and thus knowing g) up to inversion with respect to the unit circle. There are in general K -1 roots, and each one can be inverted. This gives 2 K-1 solutions.

We set

Q(g)(X) = P (g)(1/X) = ĝ[K]X -(K-1) + ĝ[K -1]X -(k-2) + ... + ĝ[1].
The equation |g| 2 = m 2 is equivalent to |g| 2 = m 2 , that is This number is the coefficient of order k of P (g)(X)Q(g)(X), so |g| = m if and only if

P (g)(X)Q(g)(X) = N K-1 k=-(K-1) m 2 [k]X k . ( 22
)
Let us denote by r 1 , ..., r K-1 the roots of P (g)(X), so that P (g)(X) = ĝ[K](X -r 1 )...(X -r K-1 );

Q(g)(X) = ĝ[K]
(1/X -r 1 )...(1/X -r K-1 ).

From ( 22), the equality |g| = m holds if and only if ĝ[K], r 1 , ..., r K-1 satisfy

|ĝ[K]| 2 K-1 j=1 (X -r j )(1/X -r j ) = N K-1 k=-(K-1) m 2 [k]X k . ( 23
)
If we denote by s 1 , 1/s 1 , ..., s K-1 , 1/s K-1 the roots of the polynomial function

K-1 k=-(K-1)
m 2 [k]X k , then the only possible choices for r 1 , ..., r K-1 are, up to permutation, r 1 = s 1 or 1/s 1 , r 2 = s 2 or 1/s 2 , . . . So there are 2 K-1 possibilities. Once the r j have been chosen, ĝ[K] is uniquely determined by ( 23), up to multiplication by a unitary complex.

From r 1 , ..., r K-1 , ĝ[K], P (g) is uniquely determined and so is g. The algorithm is summarized in Algorithm 2.

Algorithm 2 reconstruction by exhaustive search for a small problem Input: K,m 1: Compute the roots of

K-1 k=-(K-1) m 2 [k]X k .
2: Group them by pairs (s 1 , 1/s 1 ), ..., (s K-1 , 1/s K-1 ). 3: List the 2 K-1 elements (r 1 , ..., r K-1 ) of {s 1 , 1/s 1 }×...× {s K-1 , 1/s K-1 }. 4: for all the elements do 5:

Compute the corresponding ĝ[K] by (23).

6:

Compute the coefficients of P (g)(X) = ĝ[K](Xr 1 )...(X -r K-1 ).

7:

Apply an IFFT to the coefficients to obtain g. 8: end for Output: the list of 2 K-1 possible values for g.

Figure 1 .

 1 Figure 1. Example of wavelets; the figure shows ψ j for J -5 ≤ j ≤ J in the Fourier domain. Only the real part is displayed.

Figure 2 .

 2 Figure 2. ψ J , ..., ψ j+1 , ψ j (in the Fourier domain), along with ψ low j and ψ high j

Figure 3 .

 3 Figure3. Proportion of local minima, for the objective functions (6) (solid line) and (8) (dashed line), as a function of -log(r). The second proportion is independent from r. Signals were realizations of Gaussian random processes, with 500 points; the wavelets were dyadic Morlet and there was no noise.

Figure 4 .

 4 Figure 4. Four window signals, whose supports cover the interval in which errors have been detected

Figure 5 .

 5 Figure 5. For an audio signal, the reconstructed value of f ψ low 7 at the scale 7 of the multiscale algorithm, in modulus (dashed line); the solid line represents the ground true. (a) Before the error correction step (b) After the error correction step

Figure 6 .

 6 Figure 6. Mean reconstruction error (13) as a function of the noise, for an audio signal representing a human voice. (a) Maximal number of iterations per local optimization step equal to 200 (b) Maximal number equal to 2000.

Figure 7 .

 7 Figure 7. Mean reconstruction error, as a function of the size of the signal; the solid blue line corresponds to the multiscale algorithm and the dashed red one to the non-multiscale one.

Figure 8 .

 8 Figure 8. Realization of a Gaussian process (left) and modulus of its wavelet transform (right)

Figure 9 .Figure 10 .

 910 Figure 9. Mean reconstruction error as a function of the noise, for Gaussian signals of size N = 256 or 10000

Figure 11 .Figure 12 .Figure 13 .

 111213 Figure 11. Mean reconstruction error as a function of the noise, for lines extracted from images, of size N = 256 or 10000

Figure 15

 15 Figure 15. mean reconstruction error as a function of the noise, for the audio signal "Rimsky-Korsakov"

Figure 14 .Figure 16 .

 1416 Figure 14. Wavelet transforms of the audio signals (a) Rimsky-Korsakov (b) "I'm sorry"

  Figure 17. error on the signal (14) as a function of the error on the modulus of the wavelet transform (13), for several reconstruction trials; the red line y = x is here to serve as a reference (a) Gaussian signals (b) lines from images (c) sums of sinusoids (d) audio signal "I'm sorry"

Figure 19 .

 19 Figure 19. Four wavelet families. (a) Morlet (b) Morlet with dilation factor 2 1/8 (c) Laplacian (d) Gammatone

Figure 20 .

 20 Figure 20. Mean reconstruction error as a function of the noise for the four wavelet families displayed in 19. (a) Lines from images (b) Audio signal "I'm sorry"

  Figure 21. for the audio signal "I'm sorry", reconstruction error as a function of the maximal number of iterations (a) with 0.01% of noise (b) with 0.6% of noise

  be a critical point of obj 2 on the real manifold C 1 × A. There exists α ∈ C -{0}, b ∈ C such that f = αf rec + b.

  (f ) P (frec) is a constant function:P (f )(z) P (f rec )(z) = α, ∀z ∈ C, for some α ∈ C.From the definitions of P (f ) and P (f rec ),f [l] = α frec [l] ∀l = 1, . . . , [N/2].As f and f rec are analytic, we actually havef [l] = α frec [l] ∀l = 1, . . . , N -1.So the Fourier transform f -αf rec has only zero entries, except maybe in l = 0, and f -αf rec is a constant function.APPENDIX C PROOF OF LEMMA V.1 Lemma. (V.1) Let m ∈ R N and K ∈ N * be fixed. We consider the problem Find g ∈ C N s.t. |g| = m and Supp(ĝ) ⊂ {1, ..., K}.

1N

  ĝ ĝ = m 2 . For each k ∈ {-(K -1), ..., K -1}, ĝ ĝ[k] = s ĝ[k -s]ĝ[-s].

The author is with CNRS and CEREMADE, Université Paris Dauphine (e-mail: waldspurger@ceremade.dauphine.fr).

covery, maxcut and complex semidefinite programming. Mathematical Programming, 149(1-2):47-81, 2015.


Nikola Perkovic

The Electrostrong Relation

In an attempt to eliminate the Landau Pole from QED by "borrowing" asymptotic freedom from QCD, I was successful in uniting the coupling constants of the two, respectively, in the "Electrostrong relation". The Landau pole, estimated at around eV, leads to a value of

in the Electrostrong Relation, instead of infinity as was the case before. In addition, the Strong CP problem is also solved proving that there is no need for fine tuning in QCD. The first step however, was improving the measurement for the running of the strong coupling constant Q that can be tested experimentally for value of (Λ Q ) for the energy scale, QCD integration parameter Λ Q offering such a prediction for the first time.

I Introduction

The two dimensionless running constants Q and Q will be united in a single equation, the Electrostrong Relation, allowing us to both provide a solution for the Strong CP problem and remove the problem arising due to infinities in Q known as the Landau pole [START_REF] Landau | Niels Bohr and the Development of Physics[END_REF] when on a very large value of Q the fine-structure constant → ∞ therefore we will effectively make QED a mathematically complete theory. The Landau pole will be eliminated since QCD enjoys asymptotic freedom [START_REF] Gross | [END_REF], which will be "borrowed" offering us exactly what we need sinceexpansion (see Appendix 1 for a short introduction to -expansion) describes the asymptotic behavior of a denominator of convergents of continued fractions. The paper will also provide the first prediction for the value (Λ Q ) and improve the accuracy for the measurement of the running values of the strong coupling constant.

In QED vacuum polarization corrections to processes involving the exchange of virtual photons result in a Q dependence, of the effective fine structure constant Q , which is known as the running of a coupling constant which is, for the running of the fine structure constant, parameterized as:

-Q = - -∆ Q
where at the scale of the Landau pole, approximately around eV, the problem of infinities arises. This problem does not exist within QCD since in the running of the strong coupling constant Q there is asymptotic freedom.

The strong interactions do not violate CP symmetry in the manner of the weak interactions since the gauge interactions in QCD involve vector currents instead of chiral currents. Experiments too do not indicate CP violation in the QCD sector. However, if one writes the additional term in the Lagrangian density as:

ℒ = π F μ A F ̃A μ
it contains the topological term with the QCD vacuum angle , which violates CP symmetry and it can, if non-zero, lead to experimentally detectable CP violating effects. In short this is reffered to as the Strong CP problem.

II Running of 𝛂 𝐒 𝐐 and 𝛂 𝐐

The QCD Lagrangian is

ℒ Q = ∑ Ψ ̅ , (i μ ∂ μ -g s μ t A μ -m )Ψ , -F μ A F A μ (1.1) (1.2)
where μ are the Dirac -matrices, the Ψ q, are quark-field spinors for a quark of flavor f and mass m , with a color-index a that runs from a = to N = and g s = √ π is the strong gauge coupling. Let's remember that Lebesgue measure was used in the introduction toexpansion. A generalization of the Lebesgue measure for any locally compact group is known as the Haar measure [START_REF] Tornier | Haar Measures[END_REF]. We know that the simple compact Lie group we need is SU N (see the paper by G. Nagy [START_REF] Nagy | On the Haar Measure of the SU(N) Group[END_REF]). From the invariant metric ds = g dx dx with a unitarity, provided that all eigenvalues different. For matrixes W and W † we get WdW † + dWW † = and in a given set of coordinates we attain the invariant measure:

d x = √det g x ∏ dx
we parameterize the matrix W as W = exp it T therefore we conclude that:

det g ~ det Q ∏| -| >
the distribution of eigenvalues on the unitary group is given by the invariant measure:

d = ∏|exp i -exp(i )| ∏ d >
which is valid for U N groups. If we impose the constraint Σ = mod π and implement it by a δ-function, the density distribution of the eigenvalues in SU N is given by the formula above as well. For the SU group this would be of the form:

H = ∏|exp i -exp(i )| = > | | ∑ e q ∑ e q ∑ e ̅ q ∑ e q ∑ e ̅ q ∑ e ̅ q | |
We must first describe the running of Q and then proceed to connect it to the running of Q . Using -functions for QCD, a SU group theory we have for i = , , , … , n

= b π - b π - b π -⋯ (2.1) (2.2) (2.3) (2.4) (2.5)
where the QCD -function has a negative sign due to gluons carrying color charge which leads to self-interactions. At one-loop order we determine:

= -b ln

Now we introduce the -expansion for a real number q:

x = ∑ c ∞ = q -
where for all n we have c q where c doesn't have to be an integer and q is a floor function [START_REF] Hilbert | Geometry and the Imaginaion[END_REF]. We connect equations 2.6 and 2.7 which yield: where the coefficient C of the -dependent piece in the term is yet to be determined and n is the number of flavors, allows us to establish that = (∑ e q ) ∑ e q ⁄ is consistent with the Vandermonde determinant of the Haar measure H for SU [START_REF] Ovchinnkov | A Course on Lebesgue's Theory[END_REF]. We introduce v = GeV the vacuum expectation value of the Higgs field [START_REF] Buchalla | Introduction to The Standard Model[END_REF]. When a local rather than global symmetry is spontaneously broken the gauge field of the broken symmetry becomes massive, leading to what should be a Goldstone scalar to become a longitudinal mode of the massive vector, which works for both abelian and non-abelian local symmetries. In the non-abelian case, for each spontaneously broken generator T 𝑎 of the local symmetry the corresponding gauge field A 𝜇 𝑎 x becomes massive. The Higgs mechanism gives quarks mass and leaves photons and gluons massless therefore we measure (Λ Q ) for the first time, being that no other theory predicted the running value of for the Λ Q parameter:

q -= Q π + 𝒪 Λ Q Q meaning that x from
(Λ Q ) = ln ( v Λ Q ) π (2.6) (2.9) (2.7) (2.8)
This allows experimental tests to determine if the predicted value of (Λ Q ) is accurate and therefore test the -expansion method. The equation for the running value of Q is then:

Q = (Λ Q ) ln ( Q Λ Q )
where it is evident that QCD enjoys asymptotic freedom. For Q we add squared values for v, Λ Q and Q. The numerical results are provided in the table 1 bellow. All the values are in great agreement with experimental values such as provided by [START_REF][END_REF], [9], [10], [11]. Particle Data Group [12] provides the world average M z = .

. The method offered by expansion offers much desired [START_REF] Dissertori | The Determination of the Strong Coupling Constant[END_REF] accuracy and simplicity.

Using the Lyapunov exponent for the rate of exponential divergence from the initial perturbed conditions, we have ℒ x ≔ lim →∞ n -log| T ′ x | where the constant has the form = τ log -log and τ ≔ dim {x ∈ [ , : ℒ x = } where we used the Housdorff dimension [START_REF] Harris | Hausdorff dimension[END_REF]. See C. Gong, Lyapunov Exponent of Classical SU(3) Gauge Theory [START_REF] Gong | Lyapunov Exponent of Classical SU(3) Gauge Theory[END_REF] for a more detailed explanation on the relationship of the gauge coupling and the SU(3) Lyapunov exponent. From this we can deduce that in eq. (1.1) for the running of the fine-structure constant:

∆ Q = π Q
Allowing us to eliminate the Landau pole in the running of Q by establishing a link between Q and Q which eliminates the infinities since QCD enjoys asymptotic freedom by establishing:

-Q = - - π Q
The results are in good agreement with experimental values [START_REF] Bethke | The 2009 World Average of S α[END_REF], [17].

Energy [GeV] Q -Q Λ Q = . . . . . . . M z = . . . M = . . M = . . (2.10)
(2.12)

(2.11)

Table 1: Measurments for the electromagnetic and strong nuclear couplings for the values of 𝛬 𝑄 = 𝑀𝑒𝑉.

In order to prove the lack of Landau pole in my equations I provide the measurements albeit the Landau pole is not relevant for particle physics but is purely of academic interest. The value of the Landau pole is estimated roughly to be around eV, for this value Q = .

and -Q = .

where the uncertainties are high since we are dealing with such a high energy level.

When it comes to the Strong CP problem we will look at the -dependence on the vacuum energy. If we take to be the phase of the overall Lagrangian ℒ = (m u ̅u + m d ̅ d) cos ⁄ + h. c. where m is the up quark mass and is m the down quark mass and we can treat the Lagrangian term as a perturbation. The only degrees of freedom which can readily be excited in QCD are pions and as M. Dine shows in his paper TASI Lectures on The Strong CP Problem [START_REF] Dine | TASI Lectures on The Strong CP Problem[END_REF], the relation of this dependence is E = E π 𝑓 π cos where E π is the rest energy of a pion, however we will only be working with charged pions π + and π -and as such we obtain:

E π = E 𝑓 π cos = n qq ̅ Ry π Q = M π
where n qq ̅ is the number of constituent quarks/anti-quarks and for charged pions n qq ̅ ≡ and the Ry is the Rydberg unit of energy that is measured at 13.605693009(84) eV. Since the rest energy of charged pions is measured to be E π = .

MeV we can use it to measure the running of the fine structure constant -Q = M π = .

effectively solving the Strong CP problem. For other Electrostrong relations, such as for electrons, muons and their relationship with the Koide formula [START_REF] Rivero | The strange formula of Dr[END_REF] see Appendix 2. and correspond to a coupling ranging from -= -≈ to where we see an agreement for M ≈ GeV in the table where -= -= .

III Conclusion and Debate

. This means that, regardless of unification, the equations measure the running of the couplings in excellent agreement with both the experimental values and the theoretical predictions of the SM.

Appendix 1: β-expansion (2.13)

Let > be a real number and T i : [ , → [ , be the -transformation such that T i x =

x -[ x] where [x] is the largest integer that does not exceed x. Every x ∈ [ , can be uniquely expanded into a finite or infinite series:

x = x + x + x + ⋯ + x + ⋯
where x = [ x] and + x = (T i x) for all n . For any Borel set A ⊆ [ , where C is a constant that depends on such that C -A A A C where the measure of the length/energy scale is T i -invariant meaning that does not implicitly depend on and if A is Lebesgue measurable [START_REF] Ovchinnkov | A Course on Lebesgue's Theory[END_REF], then A = sup{ k : k ⊂ A} where k is compact.

|x - p x q x | > x -x
In the opposite case when log < π log ⁄

we have an integer N and for all n N :

|x - p x q x | < x -x
where the convergents of the -expansion of x are x = x ⁄ + x ⁄ + x ⁄ + ⋯ + x ⁄ .

Appendix 2: Electrostrong Relations for Electrons, Muons and the Higgs Boson

The relationship between the fine-structure constant and the electron rest energy/mass is well known and established but its relationship with the rest energy of muons, taus and the Higgs boson is not detailed enough. For the time being it has been well established that at the energy level of Z bosons, the fine-structure constant is approximately 1/128 and the trong coupling is approximately 0.1184 and its relationship with electron, muon and Higgs boson rest energy/mass is not explained at all. For the rest energy of electrons, we have an equation:

E = Q Ry s Q = Ry s
where Ry is the Rydberg unit of energy, M is the fine structure constant that equals .

, value provided by NIST, s ≡ is the strong coupling constant at zero Q momentum transfer and evidently Q ≡ is part of the Koide formula and shouldn't be confused with the momentum scale Q. For muons, we have a similar equation where Q = M μ the rest energy equals:

E μ = Q Ry s (M μ ) Q M μ = Ry s (M μ ) M μ
where we have to account for the running of the coupling constants. For the energy level of muons the fine-structure constant is approximated to be around 1/136 in which case the strong coupling constant s (M μ ) = . and evidently the second part of the Koide formula Q = . The Koide formula relates the rest energies/masses of three leptons, namely electrons, muons and taus'. It's been discovered by dr. Yoshio Koide in 1981. it states that Q = / in the formula:

E + E μ + E τ = Q Q (√E + √ E μ + √E τ ) = (√E + √ E μ + √E τ )
where the parameters Q and Q are the addition of this paper. Additionally we have a formula for rest energy of the Higgs boson that has been measured to be approximately 125.36(41) GeV [20]. As such, we have:

E 𝐻𝑖𝑔𝑔𝑠 = Ry s (M 𝐻𝑖𝑔𝑔𝑠 ) M 𝐻𝑖𝑔𝑔𝑠 (2) (3) (4)
(5)

  equation 2.7 equals x = Q which is the correction due to QCD effects and Λ Q = is the integration constant which coresponds to the scale where the perturbatively defined coupling would deverge and we use a value ofMeV but other values can be used as well. Knowing the values c = , c = .

  Standard Model, constants in the one-loop -functions are given as b = T R d(R )d R + T S d(S )d S -C G we can summarize these results and insert the number of fermion generations N = and Higgs doublets N = obtaining us b = / , b = -/ and b = -. These three intersections point towards a range from M = GeV to M = GeV

  If x ∈ [ , then any irrational number can be written asp x q x ⁄ = [c x , c x , … , c x , … ]where c x = [x -] and c + x = a T x for any n is a convergent of the continued fraction expansion of x where p x and q x are relatively prime which allows us to study the asymptotic behavior of the denominator k For any > there exists two positive constants A and which depend on thus for all n :If log > π log ⁄ then for almost all x ∈ [ , exists a positive integer N that depends on x, such that for all n N :