Laurent Condat

Atomic Norm Minimization for Decomposition into Complex Exponentials and Optimal Transport in Fourier Domain

Keywords: atomic norm, infinite dictionary, truncated moment problem, trigonometric moments, total variation norm, super-resolution, optimal transport

This paper is devoted to the decomposition of vectors into sampled complex exponentials; or, equivalently, to the information over discrete measures captured in a finite sequence of their Fourier coefficients. We study existence, uniqueness, and cardinality properties, as well as computational aspects of estimation using convex semidefinite programs. We then explore optimal transport between measures, of which only a finite sequence of Fourier coefficients is known.

Introduction

Mathematical analysis aims at understanding complex objects by expressing them in terms of elementary objects, the atoms. The discrete Fourier transform is a typical example, which allows us to express a vector of size N as a unique linear combination of N orthogonal atoms, which are the sampled complex exponentials of frequencies the multiples of 1/N . This framework is restrictive, however, and one may want to remove the constraint that the frequencies lie on this discrete grid. The modern paradigm of signal processing integrates the idea that a signal can be represented in a dictionary containing much more than N atoms, by choosing, among the multitude of possible representations offered by the redundancy, the simplest one [1]. Here, we study the representation of vectors as linear combinations of sampled complex exponentials, without any constraint on the frequencies. Thus, the dictionary contains the continuous infinity of sampled complex exponentials, with frequencies in [0, 1).

We define the unit circle T = {z ∈ C : |z| = 1} of the complex plane C. For every z ∈ C, arg(z) ∈ [0, 2π ) is the argument of z, such that z = |z|e j arg (z) , where j = √ -1. We endow T with the distance d : (z 1 , z 2 ) → min |arg(z 1 ) -arg(z 2 )|/(2π ), 1 -|arg(z 1 ) -arg(z 2 )|/(2π ) .

Let M ≥ 1 be an integer. We define the real Hilbert space

V = = ( m ) M m=-M ∈ C 2M +1 : (∀m) -m = * m , (1)
where • * indicates complex conjugation, endowed with the usual inner product

, ′ = M m=-M m ′ * m = 0 ′ 0 + 2 Re M m=1 m ′ * m ∈ R. (2)
We consider the infinite family A = a(z) z ∈T of elements of V, defined by a(z) m = z -m = e -jm arg(z) , m = -M, . . . , M.

(

We call A the dictionary and its elements a(z), which are sampled complex exponentials, the atoms [START_REF] Mallat | Matching pursuit in a time-frequency dictionary[END_REF][START_REF] Chen | Atomic decomposition by basis pursuit[END_REF][START_REF] Donoho | Optimally sparse representation in general (nonorthogonal) dictionaries via ℓ 1 minimization[END_REF][START_REF] Chandrasekaran | The convex geometry of linear inverse problems[END_REF][START_REF] Bhaskar | Atomic norm denoising with applications to line spectral estimation[END_REF][START_REF] Tang | Compressed sensing off the grid[END_REF]. The frequency of a(z) is defined as arg(z)/(2π ) ∈ [0, 1).

In the following, is an arbitrary nonzero element of V.

We focus on the decomposition of into atoms. For this, we need a parameterization of A; that is, a way to express as a linear combination of atoms. Since Acontains a non-denumerable infinity of atoms, it is natural to work with finite signed Borel measures on T, the set of which is denoted by M. So, the analysis of consists in determining a measure µ ∈ M, such that

m = ∫ T z -m dµ(z) = ∫ 1 0 e -j2π f m dµ(e j2π f ), m = -M, . . . , M. (4)
Thus, the m are trigonometric moments, or Fourier coefficients, of µ. Therefore, we write = Fµ [START_REF] Chandrasekaran | The convex geometry of linear inverse problems[END_REF] if the relation in ( 4) is satisfied, and we say that µ explains . So, one can view the problem as the identification of a measure from a truncated sequence of its trigonometric moments, the so-called truncated moment problem [START_REF] Curto | Recursiveness, positivity, and truncated moment problems[END_REF]. If one finds µ from = ( m ) M m=-M , one can calculate all its moments ( m ) m ∈Z and so extrapolate , hence the term super-resolution [START_REF] Candès | Towards a mathematical theory of superresolution[END_REF][START_REF] Azaïs | Spike detection from inaccurate samplings[END_REF][START_REF] Benedetto | Super-resolution by means of Beurling minimal extrapolation[END_REF]: one recovers high-frequency information of µ, given only low-frequency information. The two "mathematical" and "signal processing" views of the problem consist in reasoning in terms of measures with moments , or in terms of decompositions of into atoms (viewing µ as the spectrum of ), respectively; this is just a matter of defining which are the time and frequency domains. We can note that the more general setting of the recovery of a measure or a spline function from linear measurements has been considered in the literature [START_REF] Fisher | Spline solutions to L1 extremal problems in one and several variables[END_REF][START_REF] Castro | Exact reconstruction using Beurling minimal extrapolation[END_REF], with several recent contributions to this topic [START_REF] Unser | Splines are universal solutions of linear inverse problems with generalized TV regularization[END_REF][START_REF] Flinth | Exact solutions of infinite dimensional total-variation regularized problems[END_REF][START_REF] Boyer | On representer theorems and convex regularization[END_REF]. In a nutshell, regularizing an inverse problem with a convex function yields solutions that are convex combinations of a small number of atoms, which depend on the geometry of the regularizer level sets. In this paper, we focus on the specific properties of recovering a measure, when the measurements are Fourier coefficients.

In this study, we are interested in atomic decompositions of ; that is, we look at questions of existence, uniqueness, number of atoms, and more generally at the information contained in . Equivalently, we look at what can be said about a discrete measure given a finite number of its Fourier coefficients. We do not consider approximation, which consists in estimating a simple approximate atomic decomposition of , e.g. when is corrupted by noise; there is a vast literature on approximation, see [START_REF] Stoica | Spectral Analysis of Signals[END_REF][START_REF] Pereyra | Exponential data fitting[END_REF][START_REF] Condat | Cadzow denoising upgraded: A new projection method for the recovery of Dirac pulses from noisy linear measurements[END_REF][START_REF]Super-resolution of positive spikes by Toeplitz low-rank approximation[END_REF][START_REF] Andersson | A new frequency estimation method for equally and unequally spaced data[END_REF][START_REF] Duval | Exact support recovery for sparse spikes deconvolution[END_REF] and references therein.

The paper is organized as follows: In Sections 2-5, we review the properties of atomic decompositions, the close relationship between positive semidefinite Toeplitz matrices and them, and the Prony method to extract the parameters from a mixture of atoms. Many results, known but scattered in the literature or unpublished because they are straightforward, are stated as remarks. In Section 6, we establish a precise characterization of the existence and uniqueness of a measure explaining with minimal total variation norm. We explain how to determine this minimal measure numerically, by a two-step procedure: first solving a convex semidefinite program, for which we propose an iterative algorithm in Section 7; and then applying Prony's estimation method. In Section 8, we exhibit a new construction of an atomic decomposition of with 2M atoms, uniformly located on T. This construction plays an important role in Section 9, where the minimal atomic decomposition of a mixture of two atoms with opposite amplitudes is characterized. Finally, in Section 10, we define new functionals to do optimal transport between two elements of V, by doing optimal transport of the atoms composing them. Equivalently, we perform optimal transport of measures, given only a finite number of their Fourier coefficients.

Atomic decompositions and atomic norm

We will see that we can reason with discrete measures only; that is, express with a finite number of atoms. So, by anticipation, we call an atomic decomposition of , the parameters

K ∈ N\{0}, z k ∈ T all distinct, c k ∈ R\{0}, such that = K k =1 c k a(z k ). (6)
If has such an atomic decomposition, the measure µ = K k =1 c k δ z k explains , where δ z is the Dirac measure located in z ∈ T.

Remark 1. The dimension of V is 2M + 1. The rank of A, which is the largest number of linearly independent atoms, is 2M + 1 as well. Indeed, the invertibility of the Vandermonde matrices

       1 e -jω 1 • • • e -j(2M +1)ω 1 . . . . . . . . . . . . 1 e -jω 2M +1 • • • e -j(2M +1)ω 2M +1        , (7)
for every set {ω 1 , . . . , ω 2M +1 } of distinct elements of [0, 2π ), implies the invertibility of the matrices

       e jM ω 1 • • • e jω 1 1 e -jω 1 • • • e -jM ω 1 . . . . . . . . . . . . . . . . . . . . . e jM ω 2M +1 • • • e jω 2M +1 1 e -jω 2M +1 • • • e -jM ω 2M +1        , (8)
which yields the given value of the rank.

So, every family a(z 1 ), . . . , a(z 2M +1 ) of 2M + 1 distinct atoms is a basis of V, i.e. is linearly independent and spans V by linear combinations with real weights.

Remark 2. The invertibility of Vandermonde matrices also implies that the spark of A, which is the smallest number of linearly dependent atoms [START_REF] Donoho | Optimally sparse representation in general (nonorthogonal) dictionaries via ℓ 1 minimization[END_REF], is 2M + 2; A has full spark [23].

A direct consequence is the following [4, Theorem 3]: if has two different atomic decompositions, they involve at least spark(A) = 2M + 2 different atoms.

Therefore, if has an atomic decomposition with 1 ≤ K ≤ M atoms, there is no other atomic decomposition of with less or as many atoms. Moreover, we will see in Section 4 that it is easy to find this decomposition. It is in the case where does not have a decomposition with 1 ≤ K ≤ M atoms, a property which is easy to test, that it is very difficult to find the decomposition of with the smallest number of atoms (this number being at most 2M, see Remark 13).

Remark 3.

There is an infinity of atomic decompositions of with at most 2M + 1 atoms. Indeed, according to Remark 1, given 2M + 1 arbitrary distinct atoms, one can express as a linear combination of them, with real weights c k , obtained by solving a linear system. Moreover, if has a decomposition with more than 2M + 1 atoms, one can re-express all the atoms in terms of a subset of size 2M + 1 of them, to obtain a decomposition of with at most 2M + 1 atoms. Other said, there is an infinity of discrete measures, made of at most 2M + 1 Dirac measures, explaining .

Remark 4.

We can exhibit the "classical" uniform atomic decomposition

= 2M +1 k =1 c k a(e j2π k /(2M+1) ), (9)
where the sequence (c k ) 2M +1 k =1 is simply the inverse discrete Fourier transform of :

c k = 1 2M + 1 M m=-M
m e j2π km/(2M+1) , k = 1, . . . , 2M + 1. [START_REF] Azaïs | Spike detection from inaccurate samplings[END_REF] However, in most applications, is a mixture of unknown atoms, which must be identified, with frequencies living in the continuum [0, 1). Enforcing the frequencies to lie on the grid of the multiples of 1/(2M + 1), or actually on any predetermined finite set of atoms, yields undesirable effects: the reconstruction or approximation of the estimated signal is not parsimonious, with spurious oscillations; this is called spectral leakage, basis mismatch, or gridding error [START_REF] Chi | Sensitivity to basis mismatch in compressed sensing[END_REF][START_REF] Fannjiang | Coherence pattern-guided compressive sensing with unresolved grids[END_REF][START_REF] Duarte | Spectral compressive sensing[END_REF]. This phenomenon is somewhat attenuated by taking a finer and finer grid, while in doing so, with more and more coherent atoms, the estimation problems become increasingly badly conditioned and numerically unstable. This justifies to place ourselves in a continuous framework, as we do in this study, instead of in a discrete framework with a finite, even if redundant, dictionary of predetermined atoms [START_REF] Donoho | Optimally sparse representation in general (nonorthogonal) dictionaries via ℓ 1 minimization[END_REF][START_REF] Fuchs | Sparsity and uniqueness for some specific under-determined linear systems[END_REF][START_REF] Bourguignon | A sparsity-based method for the estimation of spectral lines from irregularly sampled data[END_REF]. Thus, we are interested in finding, among the infinite multitude of possible atomic decompositions of , the "simplest" ones, with respect to some notion of parsimony. Finding the atomic decomposition with the smallest number of atoms is a very difficult nonconvex problem. That is why we focus on its convex relaxation and consider the form of optimality defined in the following.

We define the atomic norm of as

a := inf µ TV : µ ∈ M, = Fµ , (11)
Strictly speaking, we should remove from this sum the indices k, for which c k = 0.

where • TV is the total variation norm of a measure (not to be confused with the total variation of an image [START_REF] Condat | Discrete total variation: New definition and minimization[END_REF]). In that respect, let us recall the following property [START_REF] Fisher | Existence, uniqueness, and minimality of the Jordan measure decomposition[END_REF]: every measure µ ∈ M has a unique Jordan decomposition (µ + , µ -) of two finite positive measures µ + and µ -, such that µ = µ +µ -, satisfying the following optimality property: for every pair of positive measures (ν + , ν -) with µ = ν +ν -, then ν + ≥ µ + and ν -≥ µ -. Hence, we can define µ TV = µ + (T) + µ -(T).

As easily checked, the atomic norm is indeed a norm on V.

We call a minimal atomic decomposition of , an atomic decomposition of , such that

K k =1 |c k | =
a . Remark 5. The atomic norm is defined in several papers [START_REF] Chandrasekaran | The convex geometry of linear inverse problems[END_REF][START_REF] Bhaskar | Atomic norm denoising with applications to line spectral estimation[END_REF][START_REF] Tang | Compressed sensing off the grid[END_REF] as a = inf t ≥ 0 : ∈ t conv(A ∪ -A) , [START_REF] Fisher | Spline solutions to L1 extremal problems in one and several variables[END_REF] where conv denotes the convex hull, which, since V is of dimension 2M + 1, is the set of convex combinations of 2M + 1 signed atoms. Therefore, we can equivalently write [START_REF] Fisher | Spline solutions to L1 extremal problems in one and several variables[END_REF]

as a = inf K k =1 |c k | : = K k =1 c k a(z k ), for some K ≤ 2M + 1, distinct z k ∈ T, c k ∈ R\{0} . (13)
We will see in Proposition 2 that these two definitions [START_REF] Benedetto | Super-resolution by means of Beurling minimal extrapolation[END_REF] and ( 12)-( 13) coincide, and that the infimum in ( 11), ( 12), ( 13) is attained; that is, it is a minimum. In other words, a minimal atomic decomposition of exists.

Relation to Toeplitz matrices

In order to characterize more precisely the atomic norm and the atomic decompositions of minimizing it, we introduce the linear operator T, which maps to the Hermitian Toeplitz matrix

T( ) :=          0 1 • • • M -1 0 • • • M-1 . . . . . . . . . . . . -M -M+1 • • • 0          . ( 14
)
We can note that T defines a bijection between V and the set of (M + 1) × (M + 1) Hermitian Toeplitz matrices, thereafter denoted by T(V). So, we will talk equally about an atomic decomposition of or of T( ).

T( ) is said to be indefinite if it is neither positive semidefinite nor negative semidefinite; that is, if it has at least one positive and one negative eigenvalues. T( ) is said to be singular if it is not invertible; that is, if it has at least one eigenvalue equal to zero.

We remark that, for every z ∈ T,

T(a(z)) =          1 z -1 • • • z -M z 1 • • • z -M+1 . . . . . . . . . . . . z M z M-1 • • • 1          =          1 z . . . z M          [1 z -1 • • • z -M ] (15)
is positive semidefinite and of rank 1.

Remark 6. If has an atomic decomposition

= K k =1 c k a(z k ) with 1 ≤ K ≤ M + 1 atoms, the rank of T( ) is K. Indeed, T( ) =          1 • • • 1 z 1 • • • z K . . . . . . . . . z M 1 • • • z M K          ×        c 1 0 . . . 0 c K        ×        1 z -1 1 • • • z -M 1 . . . . . . . . . . . . 1 z -1 K • • • z -M K        (16)
and the rank of A H BA, where A is a matrix with full column rank and • H denotes the Hermitian transpose, is equal to the rank of B [31, Lemma 1].

As we see, the rank of T( ) will play an important role in the analysis. At this stage, one can wonder if the converse of the last property is true; that is, if there always exists an atomic decomposition of with K = rank(T( )) atoms. Unfortunately, this is not the case, as we will see in Section 4.

Remark 7.

One can think that if has an atomic decomposition with M + 1 < K ≤ 2M + 1 atoms, the rank of T( ) is M + 1. This is not true in general and we can construct examples where the rank is as low as

2M + 2 -K. Indeed, if 2M +2 k =1 c k a(z k ) = 0 (given distinct z k , such c k 0 exist, since the a(z k ) are linearly dependent), then for every M + 1 < K ≤ 2M + 1, = K k =1 c k a(z k ) = 2M +2 k =K +1 (-c k )a(z k )
, and there are two different atomic decompositions of , the second one having less atoms. The rank of T( ) is then 2M + 2 -K.

Remark 8. The trace of T( ), which is also the sum of its eigenvalues, is equal to (M + 1) 0 . So, it is equal to (M + 1)µ(T), for every measure µ explaining , and to (M + 1) K k =1 c k , for every atomic decomposition of . As we see, decomposition into exponentials is intimately related to the spectral properties of Hermitian Toeplitz matrices.

Prony's estimation method

The following estimation method, already known by G. R. de Prony in the 18th century [START_REF] Riche De Prony | Essai expérimental et analytique : sur les lois de la dilatabilité de fluides élastiques et sur celles de la force expansive de la vapeur de l'eau et de la vapeur de l'alcool à différentes températures[END_REF], allows us to know, in a simple way, whether has an atomic decomposition with K ≤ M atoms, and if so, to find the corresponding parameters. We outline the method in this section, see [START_REF] Auton | Investigation of procedures for automatic resonance extraction from noisy transient electromagnetics data[END_REF][START_REF] Blu | Sparse sampling of signal innovations[END_REF][START_REF] Potts | Nonlinear approximation by sums of nonincreasing exponentials[END_REF] for more details.

First, let us suppose that there exists a decomposition

= K k =1 c k a(z k ) with 1 ≤ K ≤ M. We define the annihilating polynomial H(z) = K k =0 h k z k := K k =1 (z -z k ).
This polynomial owes its name to the fact that it cancels the sequence by convolution [START_REF] Blu | The generalized annihilation property-A tool for solving finite rate of innovation problems[END_REF][START_REF] Dragotti | Sampling moments and reconstructing signals of finite rate of innovation: Shannon meets Strang-Fix[END_REF]:

K k =0 h k m-k = K k =0 h k K l =1 c l z k -m l = K l =1 c l z -m l H(z l ) = 0, m = -M + K, . . . , M. (17)
Let us form the rectangular Toeplitz matrix of size (2M

+ 1 -K) × (K + 1), R K ( ) :=          M-K • • • M M-K -1 • • • M-1 . . . . . . . . . -M • • • -M+K          . ( 18
)
Then, the equations ( 17) can be rewritten as

R K ( )[h K • • • h 0 ] T = [0 • • • 0] T . (19)
Thus, the vector

[h K • • • h 0 ] T is in the kernel of R K ( ). More precisely, R K ( ) is of rank K and its kernel, of dimension 1, is spanned by [h K • • • h 0 ] T .
Now, we consider the converse of these properties. Let K be the rank of T( ), with K ≤ M (if T( ) is positive definite, see Remark 11, whereas if T( ) is of full rank and is not positive definite, the method cannot be applied to find an atomic decomposition of ). Let us form the matrix R K ( ) as in [START_REF] Pereyra | Exponential data fitting[END_REF]. It is of rank K as well. Let [h K • • • h 0 ] T be the unique, up to a constant, nonzero vector in the kernel of R K ( ). For a numerically robust estimation, we compute the singular value decomposition of R K ( ) and we take [h K • • • h 0 ] T as the right singular vector corresponding to the smallest singular value. Then, we compute the roots z k of the polynomial H(z

) := K k =0 h k z k . If we find K distinct roots belonging to T, then has an atomic decomposi- tion = K k =1 c k a(z k )
, where the amplitudes c k are obtained by solving the least-squares linear system

Q H Q[c 1 • • • c K ] T = Q H [ -M • • • M ] T , (20)
where

Q :=        z -M 1 • • • z -M K . . . . . . . . . z M 1 • • • z M K        . ( 21
)
In the other cases, does not have an atomic decomposition with at most M atoms. Note that the annihilating polynomial can have less than K nonzero roots. For instance, for

= (2, 1, 1, 1, 1, 1, 2), T( ) =         1 1 1 2 1 1 1 1 1 1 1 1 2 1 1 1         (22)
is of rank K = 3, so R K ( ) = T( ) and the annihilating polynomial is H(z) = z(z -1), which has only one nonzero root.

In the case where H(z) has K distinct nonzero roots z k , but not all in T, we can write

m = K k =1 c k z -m
k , but this amounts to decomposing into other objects than the atoms of A. This is linked to the fact that the Hermitian Toeplitz matrices of rank 1 are of the form T((z -M , . . . , z M )), but for an arbitrary z ∈ C\{0}, not necessarily in T.

We will see in the following that if T( ) is positive semidefinite, Prony's method always yields K distinct roots in T. If T( ) is indefinite, the key for decomposing will consist in splitting T( ) into the sum of two positive semidefinite and negative semidefinite Toeplitz matrices.

Case of positive decompositions

For every measure µ explaining , µ TV ≥ µ + (T) ≥ µ + (T)µ -(T) = 0 , so we always have a ≥ 0 . Moreover, if a positive measure µ explains , we have µ TV = µ(T) = 0 , so that the infimum in [START_REF] Benedetto | Super-resolution by means of Beurling minimal extrapolation[END_REF] is attained for this measure, with a = µ TV = 0 . So, it is interesting to address, in the first place, the case where is explained by a positive measure.

Remark 10. If a discrete measure K k =1 c k δ z k explains , with 1 ≤ K ≤ M + 1 and c k > 0,
then T( ) is positive semidefinite and of rank K. Indeed, ∀z ∈ T, T(a(z)) is positive semidefinite and the set of positive semidefinite matrices is a closed convex cone. More precisely, for every

vector x ∈ C M+1 , x H T( )x = ∫ T x H T(a(z))x dµ(z) = ∫ T | M+1 m=1 z 1-m x m | 2 dµ(z) ≥ 0.
Moreover, we have seen in Remark 6 that T( ) is of rank exactly K.

A major result, which implies the converse of the property in Remark 10, dates back to the work of Carathéodory [START_REF] Carathéodory | Über den Variabilitätsbereich der Fourierschen Konstanten von positiven harmonischen Funktionen[END_REF][START_REF] Carathéodory | Über den Zusammenhang der Extremen von harmonischen Funktionen mit ihren Koeffizienten und über den Picard-Landau'schen Satz[END_REF]. It can be formulated as follows:

Carathéodory's Theorem. T( ) is positive and of rank K, with 1 ≤ K ≤ M + 1, if and only if there exist distinct elements z k ∈ T and positive reals c k , such that

= K k =1 c k a(z k ). Moreover, if K ≤ M, the discrete measure K k =1 c k δ z k is the unique positive measure explaining . Thus, if T( ) is positive semidefinite, a positive discrete measure µ = K k =1 c k δ z k , with K ≤ M +1, explains , and we have µ TV = K k =1 c k = 0 = a .
In other words, the Theorem states the existence of a minimal atomic decomposition of , if T( ) is positive semidefinite.

Moreover, according to Remark 2, there exists no decomposition of with less atoms, so the obtained atomic decomposition is doubly optimal.

Remark 11. If T( ) is positive definite, with rank K = M + 1, there exists an infinity of minimal atomic decompositions of with M + 1 atoms. We can construct all of them as follows [START_REF] Yang | Exact joint sparse frequency recovery via optimization methods[END_REF]Remark 2.1]: we choose z M+1 ∈ T and we set

c M+1 := [1 z -1 M+1 • • • z -M M+1 ]T( ) -1 [1 z M+1 • • • z M M+1 ] T -1 . ( 23
)
Then c M+1 > 0 and

T -c M+1 a(z M+1 ) = T( ) -[1 z M+1 • • • z M M+1 ] T [1 z -1 M+1 • • • z -M M+1 ] (24)
is positive semidefinite of rank M. Carathéodory's Theorem applied to this residual yields -

c M+1 a(z M+1 ) = M k =1 c k a(z k ).
All in all, we end up with the positive minimal atomic decomposition = M+1 k =1 c k a(z k ). For instance, for = (0, . . . , 0, 1, 0, . . . , 0), i.e. T( ) = Id, the identity matrix, we get the atomic decompositions, for every φ ∈ [0, 2π ),

= 1 M + 1 M+1 k =1
a e j(2π k +φ)/(M+1) .

(25)

Remark 12. There is a complete proof of Carathéodory's Theorem by Curto and Fialkow [8, section 6]. Given Remarks 10 and 11, there only remains to prove that, if T( ) is positive semidefinite of rank 1 ≤ K ≤ M, then a positive measure explaining exists, is unique, and is concentrated at K points of T. The difficult part, which is not detailed here, is to show the existence. The rest of the proof is simple: we consider the annihilating polynomial

H(z) = K k =0 h k z k = K k =1 (z -z k ) defined in Section 4 and we define the vector of size M + 1, h = [h K • • • h 0 0 • • • 0] T . Let µ ∈ M be a positive measure explaining . We have 0 = hH T( ) h = 0≤k ,l ≤K h * k h l k -l = ∫ T 0≤k ,l ≤K h * k h l z l -k dµ(z) (26) = ∫ T |H(z)| 2 dµ(z).
Therefore, µ is a discrete measure concentrated at the roots z k of H(z) belonging to T. Since µ cannot be made of less than K Diracs measures (else the rank of T( ) would be less than K, according to Remark 10), all the z k , k = 1, . . . , K, are in T. Finally, the unique coefficients c k such that = K k =1 c k a(z k ), defined in [START_REF]Super-resolution of positive spikes by Toeplitz low-rank approximation[END_REF], are positive, according to Remark 9. Corollary 1. If a positive measure, which is not a discrete measure concentrated in at most M points, explains , then T( ) is of rank M + 1. This is specific to the positive case, see Remark 7 for the general case.

Corollary 2. An atom cannot be expressed as a linear combination of other atoms with positive coefficients. So, the atoms are the extremal points of the convex hull of A.

General case

Remark 13. In the previous section, we have seen that if T( ) is positive semidefinite, has a decomposition with at most M + 1 atoms. The analysis is of course the same when T( ) is negative semidefinite. So, let us direct our attention to the case where T( ) is indefinite. We have seen in Remark 3 that there are atomic decompositions of with 2M + 1 atoms. We can now improve this bound and construct decompositions of with at most 2M atoms. We will see in Proposition 3 one such decomposition. Another construction is the following. Let i -< 0 be the smallest eigenvalue of T( ). T( )i -Id is positive semidefinite and singular, so by Carathéodory's Theorem, it can be decomposed into at most M atoms. As in Remark 11, we can decompose i -Id into M + 1 atoms, one of which being an atom of T( )i -Id. By combining them, T( ) = (T( )i -Id) + i -Id can indeed be decomposed into at most 2M atoms. Remark 14. The bound above of 2M atoms is optimal: it is possible that has a decomposition with 2M atoms and not less. We can give the example [41, Example 3.4], for M = 2, of

= (3, 1, 1, 1, 3) = 3 2 a(1) - 1 2 a(j) + 1 2 a(-1) - 1 2 a(-j), (27)
for which it is shown that no decomposition with less than four atoms exists [START_REF] Bakonyi | Moment problems for real measures on the unit circle[END_REF].

We now focus on minimal atomic decompositions of . For this purpose, let us consider the three equivalent formulations of the same convex optimization problem, where tr denotes the trace, 0 denotes positive semidefiniteness, and V := T( ) :

V + := argmin X tr(X ) s.t. X ∈ T(V) and X 0 and X -V 0, (28)

≡ V + := argmin X 1 M+1 tr(X ) + tr(X -V ) (29)
s.t. X ∈ T(V) and X 0 and X -V 0, ≡ V + := T( + ) , where + := argmin

x ∈V (2x 0 -0 ) s.t. T(x) 0 and T(x -) 0. ( 30
)
Proposition 1. The solution V + of the problem (28) exists and is unique. Accordingly, let us define

V -:= V + -V . Moreover, if V is positive semidefinite, V + = V and V -= 0; if V is negative semidefinite, V + = 0 and V -= -V ; else, V + and V -are singular.
Proof: the case where V is positive semidefinite or negative semidefinite is obvious, so let us suppose that V is indefinite. We first remark that the feasible set of the problem ( 28) is not empty: let i -< 0 be the smallest eigenvalue of V ; X = Vi -Id is positive semidefinite, Hermitian and Toeplitz, as well as X -V .

Then the trace defines a norm for the positive semidefinite matrices, and, as a consequence of the Weierstrass theorem [42, Theorem 1.28], a norm attains its minimum on a closed set in finite dimension, so a solution V + exists.

Let V + be a solution of [START_REF] Bourguignon | A sparsity-based method for the estimation of spectral lines from irregularly sampled data[END_REF]. Set V -:= V + -V . Suppose that V + is positive definite. Then there exists 0 < a < 1, such that X := V + -aV -is positive semidefinite. Moreover, X ∈ T(V) and X -V = (1a)V -is positive semidefinite. But tr(X ) < tr(V + ), which contradicts the optimality of V + . Hence, V + is singular. By the same reasoning, V -is singular as well.

Let us finally show that the solution of ( 28) is unique. Let X 1 and X 2 be two solutions. Then X = (X 1 + X 2 )/2 is a solution as well, so X is singular. By Carathéodory's Theorem, X has an atomic decomposition

K k =1 c k a(z k ), where 1 ≤ K ≤ M is the rank of X . Let us define the annihilating polynomial H(z) = K k =0 h k z k := K k =1 (z -z k ) and the vector of size M + 1, h = [h K • • • h 0 0 • • • 0] T .
We have hH X h = 0, so hH X 1 h + hH X 2 h = 0, and since X 1 and X 2 are positive semidefinite, hH X 1 h = hH X 2 h = 0. As in [START_REF] Duarte | Spectral compressive sensing[END_REF], this implies that the positive atomic decompositions of X 1 and X 2 are of the form

K k =1 c 1,k a(z k ) and K k =1 c 2,k a(z k )
, respectively, for some c 1,k ≥ 0 and c 2,k ≥ 0. By the same reasoning, X 1 -V and X 2 -V have atomic decompositions of same support, of the form

K ′ k =1 c ′ 1,k a(z ′ k ) and K ′ k =1 c ′ 2,k a(z ′ k ), respectively, with 1 ≤ K ′ ≤ M, some z ′ k ∈ T all distinct and different from the z k , and some c ′ 1,k ≥ 0 and c ′ 2,k ≥ 0. Therefore, = K k =1 c 1,k a(z k ) -K ′ k =1 c ′ 1,k a(z ′ k ) = K k =1 c 2,k a(z k ) -K ′ k =1 c ′ 2,k a(z ′ k ). By subtraction, K k =1 (c 1,k -c 2,k )a(z k ) + K ′ k =1 (c ′ 2,k -c ′ 1,k )a(z ′ k ) = 0 V
, the zero sequence of V. Since the atoms are linearly independent, for every k, c 1,k = c 2,k and c ′ 1,k = c ′ 2,k . Hence, X 1 = X 2 and the solution, call it V + , is unique.

Consequently, to every V ∈ T(V), we can assign a unique pair (V + , V -) of Hermitian Toeplitz positive semidefinite matrices, with V + solution of (28) and V -= V + -V . By analogy with the Jordan decomposition of measures, we call (V + , V -) the Jordan decomposition of V . Equivalently, to every ∈ V, we can assign a unique pair ( + , -) ∈ V 2 , with + solution of (30) and -= + -. We call ( + , -) the Jordan decomposition of .

; that is + = K + k =1 c k a(z k ) and -= K + +K - k =K + +1 (-c k )a(z k )
, where K = K + + K -, K + := rank(V + ), K -:= rank(V -), V := T( ). The existence of these two atomic decompositions is guaranteed by Carathéodory's Theorem, since V + 0 and V -0. Thus,

a = µ ⋆ TV = K k =1 |c k | = 1 M+1 tr(V + ) + tr(V -) = + 0 + - 0 = 2 + 0 -0 . (32)
µ ⋆ = K k =1 c k δ z k exhibited in (i)
is the unique solution of the problem [START_REF] Satorra | A theorem on the rank of a product of matrices with illustration of its use in goodness of fit testing[END_REF]. Then we have 0 ≤ K + ≤ M and 0 ≤ K -≤ M. On the contrary, if T( ) is positive definite (resp. negative definite), there is an infinity of measures, necessarily positive (resp. negative), solution of [START_REF] Satorra | A theorem on the rank of a product of matrices with illustration of its use in goodness of fit testing[END_REF], among which an infinity of discrete measures concentrated at exactly K = M + 1 points of T; we can construct them explicitly, according to Remark 11. Proof: (i) let µ ∈ M, with Jordan decomposition (µ + , µ -), explaining . Let (x + , x -) ∈ V 2 , such that µ + explains x + and µ -explains x -. We have x +x -= and µ TV = x + 0 + x - 0 = 2x + 0 -0 . Since x + is feasible for the problem (30), the solution + of which exists and is unique, µ is a solution of (31) if and only if x + = + (and x -= -). The rest of the statement follows. (ii) Suppose that V = T( ) is neither positive definite nor negative definite. Then V + and V -are singular, according to Proposition 1, so K + = rank(V + ) ≤ M and K -= rank(V -) ≤ M. Let µ be a solution of [START_REF] Fisher | Existence, uniqueness, and minimality of the Jordan measure decomposition[END_REF]. µ + and µ -explain + and -, respectively, as shown above. But according to Carathéodory's Theorem, such positive measures are unique. So, the construction in (i) yields the unique solution.

So, we have shown that the infimum is attained in ( 11), ( 12), [START_REF] Castro | Exact reconstruction using Beurling minimal extrapolation[END_REF]. Thus, it is licit to replace the infimum by a minimum in these definitions of the atomic norm. In addition, we can strengthen the constraint K ≤ 2M + 1 to K ≤ 2M in [START_REF] Castro | Exact reconstruction using Beurling minimal extrapolation[END_REF].

Finally, we can write, for every ∈ V,

a = min X 2
M+1 tr(X ) -0 s.t. X ∈ T(V) and X 0 and X -T( ) 0. [START_REF] Auton | Investigation of procedures for automatic resonance extraction from noisy transient electromagnetics data[END_REF] We can note that this semidefinite characterization of the atomic norm was given in [43, eq. 12], without any proof.

If we go back to the problem of finding a decomposition of with the smallest number of atoms, it amounts to: minimize

X rank(X ) + rank(X -V ) s.t. X ∈ T(V)
and X 0 and X -V 0.

(34) Such nonconvex problems of rank minimization are notoriously difficult, and atomic norm minimization is a convex and efficient alternative.

Given an atomic decomposition

= K k =1 c k a(z k )
, we define its separation D({z k }) := min k l d(z k , z l ), as the minimal distance between its atoms. A remarkable result by Candès and Fernandez-Granda [START_REF] Candès | Towards a mathematical theory of superresolution[END_REF] is the following: if has an atomic decomposition with large enough separation, then it is both the minimal atomic decomposition of and the atomic decomposition of with the fewest atoms; that is, the solutions of ( 29) and ( 34) are the same and correspond to this decomposition. A sufficient condition for this is D({z k }) ≥ 1.87/M (and M ≥ 128) [START_REF] Candès | Towards a mathematical theory of superresolution[END_REF]. This condition was further strengthened to D({z k }) ≥ 1.26/M (and M ≥ 1000) [START_REF] Fernandez-Granda | Super-resolution of point sources via convex programming[END_REF]. Some authors conjecture that D({z k }) ≥ 1/M is a sufficient condition [START_REF] Duval | Exact support recovery for sparse spikes deconvolution[END_REF], see also [START_REF] Da | A tight converse to the spectral resolution limit via convex programming[END_REF]. Note that every minimal atomic decomposition has a distance d(z k , z l ) ≥ 1/(2M) between two atoms, if c k is positive and c l is negative [22, Corollary 2]. Of course, these sufficient conditions imply that the number of atoms is K ≤ M, and we have seen that in that case, Prony's method yields this optimal decomposition directly. But we cannot hope for a stable method to recover atomic decompositions with K ≥ M + 1 atoms, characterized by 2K ≥ 2M + 2 real parameters, from , characterized by 2M + 1 real degrees of freedom. More importantly, the fact that separated enough atomic decompositions are represented by solutions of convex optimization problems, opens the door to their robust approximate estimation by convex minimization, from inexact or noisy coefficients m ≈ m . For instance, given ∈ V, we can minimize

x ∈V -x 2 2 + λ x a , (35)
≡ minimize

(x + ,x -)∈V 2 -x + + x -2 2 + λ(x + 0 + x - 0 ) s.t. T(x + ) 0 and T(x -) 0, (36)
for some regularization parameter λ > 0. Approximation goes beyond the scope of this paper and we refer to [START_REF] Hofmann | A convergence rates result for Tikhonov regularization in Banach spaces with non-smooth operators[END_REF] and [START_REF] Bredies | Inverse problems in spaces of measures[END_REF] for some results on the topic.

Remark 15.

A measure µ ∈ Mis said to be identifiable if, given its sequence of Fourier coefficients = ( m ) M m=-M , it is the unique solution of the problem [START_REF] Satorra | A theorem on the rank of a product of matrices with illustration of its use in goodness of fit testing[END_REF]. That is, an identifiable measure can be perfectly recovered from = Fµ. 7 Numerical resolution of [START_REF] Fisher | Existence, uniqueness, and minimality of the Jordan measure decomposition[END_REF] Since the equivalent problems ( 28) and ( 30) are semidefinite programs, many efficient algorithmic strategies exist to solve them to an arbitrary precision. It is beyond the scope of this paper to compare them. We present a first order iterative algorithm, which is simple to implement and turns out to be fast. It is an instance of the overrelaxed version [START_REF] Condat | A primal-dual splitting method for convex optimization involving Lipschitzian, proximable and linear composite terms[END_REF][START_REF] Condat | Proximal splitting algorithms: Relax them all![END_REF] of the Chambolle-Pock algorithm [START_REF] Chambolle | A first-order primal-dual algorithm for convex problems with applications to imaging[END_REF]. The algorithm is the following, where T * is the adjoint operator of T and P 0 denotes the projection onto the convex cone of negative semidefinite matrices, by computing the eigendecomposition and setting the positive eigenvalues to zero.

Algorithm 1. Input: ∈ V. Output: estimate x (i) of + , the solution of [START_REF] Fisher | Existence, uniqueness, and minimality of the Jordan measure decomposition[END_REF]. Choose the proximal parameter τ > 0 and the relaxation parameter ρ ∈ [1, 2). Choose the initial estimate

x (0) ∈ V of + . Set σ := 1/τ /(2 T 2 ) = 1/τ /(2M + 2). Set U (0) 1 := U (0) 2 := T(0 V ), the zero matrix. Set V := T( ). Then iterate, for i = 0, 1, . . . 1. x(i+1) := x (i) -τ T * (U (i) 1 + U (i) 2 ) -τ (1, 0, . . . , 0), 2.
U (i+1)

1 := P 0 U (i) 1 + σ T(2x (i+1) -x (i) ) , 3. U (i+1) 2 := P 0 U (i) 2 -σV + σ T(2x (i+1) -x (i) ) , 4. (x (i+1) , U (i+1) 1 , U (i+1)
2

) := ρ(x (i+1) , U (i+1) 1 , U (i+1) 2
)

+ (1 -ρ)(x (i) , U (i) 1 , U (i) 2 ).
The algorithm converges to the solution + of ( 30) [START_REF] Condat | A primal-dual splitting method for convex optimization involving Lipschitzian, proximable and linear composite terms[END_REF]Theorem 5.3]. Moreover, the number of positive eigenvalues of

U (i) 1 + σ T(2x (i+1) -x (i) ) and U (i) 2 -σV + σ T(2x (i+1) -x (i)
) are robust estimates of K + = rank(V + ) and K -= rank(V -), the number of positive and negative atoms in the minimal atomic decomposition of , respectively. The complexity of every iteration is dominated by the computation of two eigendecompositions, so the algorithm has complexity O(M 3 ) per iteration. In practice, we initialize the algorithm with the positive part of the decomposition given in Proposition 3; that is, given the parameters {z k }, {c k } of this decomposition, we set x (0) = 2M k =1 max(c k , 0)a(z k ). So, x (0) is feasible for the problem [START_REF] Fisher | Existence, uniqueness, and minimality of the Jordan measure decomposition[END_REF]; that is, T(x (0) ) 0 and T(x (0) -) 0.

We can note that the overrelaxed Chambolle-Pock algorithm can also be used to solve [START_REF] Blu | The generalized annihilation property-A tool for solving finite rate of innovation problems[END_REF].

An alternative way to compute the atomic norm of is to solve the semidefinite program [START_REF] Bhaskar | Atomic norm denoising with applications to line spectral estimation[END_REF][START_REF] Tang | Compressed sensing off the grid[END_REF][START_REF] Candès | Towards a mathematical theory of superresolution[END_REF] minimize

q ∈C 4M +1 s.t. q -m =q * m , t ∈R 1 2 (q 0 + t) s.t.          -M T(q)
. . .

M M • • • -M t          0, ( 37
)
the optimal value of which is equal to a . Since the constraint in (37) implies T(q) 0, Carathéodory's Theorem is applicable and Prony's method applied to T(q) allows us to recover the parameters z k of the minimal atomic decomposition of . It is easy to show that at optimality, t = q 0 = a , so that it is better, for faster convergence, to remove the unnecessary variable t and to consider the problem minimize

q ∈C 4M +1 s.t. q -m =q * m q 0 s.t.          -M T(q) . . . M M • • • -M q 0          0. (38)
A typical algorithm to solve [START_REF] Dragotti | Sampling moments and reconstructing signals of finite rate of innovation: Shannon meets Strang-Fix[END_REF] or [START_REF] Carathéodory | Über den Variabilitätsbereich der Fourierschen Konstanten von positiven harmonischen Funktionen[END_REF] needs one eigendecomposition per iteration, instead of two with the proposed Algorithm 1, but it involves a matrix twice larger. Since the complexity of the eigendecomposition scales with the cube of the matrix size, we expect Algorithm 1 to be faster. We leave for future work a numerical comparison of different formulations and algorithms to compute the atomic norm of .

Uniform decomposition with at most 2M atoms

In this section, we present a new construction of a uniform atomic decomposition of , with at most 2M atoms.

Proposition 3.

There exists an atomic decomposition of , with at most 2M atoms located uniformly on T, with an explicit form:

= 2M k =1 c k a(z k ), (39)
where z k := e j(2π k /(2M)-ϕ ) and ϕ := arg( M )/M

(if M = 0, every ϕ is suitable), and the c k are obtained by inverse discrete Fourier transform from the sequence ( m e -jmϕ ) M-1 m=-M :

c k := 1 2M M-1 m=-M ( m e -jmϕ )e j2π mk /(2M) , k = 1, . . . , 2M. (41)
Proof: we first remark that the c k in (41) are real, since -m e -j(-m)ϕ = ( m e -jmϕ ) * and this value is real for m = 0 and m = -M. Then it is sufficient to notice that

m e -jmϕ = 2M k =1 c k e -j2π mk /(2M) , m = -M, . . . , M -1. (42)
So (we check that the formula is also valid for m = M),

m = 2M k =1 c k e -jm(2π k /(2M)-ϕ ) , m = -M, . . . , M. (43)
Note that the number of atoms in the proposition is at most 2M, and not exactly 2M, because some c k can be zero in [START_REF] Carathéodory | Über den Zusammenhang der Extremen von harmonischen Funktionen mit ihren Koeffizienten und über den Picard-Landau'schen Satz[END_REF] and [START_REF] Bakonyi | Moment problems for real measures on the unit circle[END_REF].

In terms of extrapolation of

= ( m ) M m=-M to the infinite sequence ( m ) m ∈Z of Fourier coefficients of µ = 2M k =1 c k δ e j(2π k /(2M )-ϕ) , we get m = e j(m-[m])ϕ [m] , m ∈ Z, (44)
where [m] = rem(m + M, 2M) -M and rem(n, N ) ∈ {0, . . . , N -1} is the remainder of the Euclidean division of n ∈ Z by N .

Remark 16. Duval and Peyré showed that there exists a measure with separation 1/(2M), which is not identifiable, by exhibiting δ e -j2π /(2M ) + δ e 0δ e j2π /(2M ) and proving that it is not the solution of ( 31) [START_REF] Duval | Exact support recovery for sparse spikes deconvolution[END_REF]Corollary 1]. According to Proposition 3, measures with separation 1/(2M), which are not identifiable, are the rule rather than the exception. Indeed, every ∈ V is explained by the measure with separation 1/(2M) constructed in Proposition 3, and there is no reason that this measure has minimal total variation norm, in general. Now, let us consider the case where is made of a single atom: = ca(z), for some c > 0 and z ∈ T (the case c < 0 is addressed the same way). T( ) is positive definite, of rank 1, so the analysis of Section 5 shows that = ca(z) is the unique positive atomic decomposition of and the unique minimal atomic decomposition of . Moreover, this single-atom decomposition coincides with the one of Proposition 3: the c k in ( 41) are all zero, except c k 0 = c for k 0 = {2M if ω = 0, 2M(ϕ + ω)/(2π ) else}, where ω = arg(z). Indeed, if ω > 0, ϕ = arg(e -jM ω )/M, so ϕ = -ω + 2πk 0 /(2M) for some even k 0 ∈ {2, 4, . . . , 2M}. Since the sequence ( m e -jmϕ = ce -j2π mk 0 /(2M) ) M-1 m=-M is a sampled complex exponential of frequency multiple of 1/(2M), its inverse discrete Fourier transforms gives c k all equal to zero, except c k 0 = c.

After this simple case of a single atom, we direct our attention to the case of two atoms with opposite amplitudes, in the next section.

Case of 2 atoms with opposite amplitudes

Let us now consider that is the sum of two atoms with opposite amplitudes. So, let z 1 and z 2 be two distinct elements of T and let c be a positive real. We set

= ca(z 1 ) -ca(z 2 ). (45)
Set ω 1 := arg(z 1 ), ω 2 := arg(z 2 ). Then we have

m = ce -jmω 1 -ce -jmω 2 , m = -M, . . . , M. (46)
We introduce the 2π -periodic function

ϒ : ω ∈ R → sin ω(2M + 1)/2 + sin ω(2M -1)/2 4M sin(ω/2) (47) = 1 4M e -jM ω + 1 2M M-1 m=-M+1 e jmω + 1 4M e jM ω , (48)
Figure 1: Plot of the function ϒ(ω), defined in [START_REF] Bredies | Inverse problems in spaces of measures[END_REF], for different values of M.

extended by continuity to ϒ(0) = 1. Note that ϒ(2πk/(2M)) = 0 for every k = 1, . . . , 2M -1. The function ϒ is plotted in Figure 1 for the first values of M.

Proposition 4. In the present case of two atoms with opposite amplitudes, the decomposition of Proposition 3 can be rewritten as (up to a circular permutation on the indices):

= 2M k =1 c k a(e jω k ), (49)
where

ω k = 2πk 2M - 2π 4M + ω 1 + ω 2 2 , (50)
c k = cϒ 2πk 2M - 2π 4M + ω 2 -ω 1 2 -cϒ 2πk 2M - 2π 4M + ω 1 -ω 2 2 . (51)
Proof: By combining ( 46), [START_REF] Bakonyi | Moment problems for real measures on the unit circle[END_REF] and the fact that -M e jMϕ = M e -jMϕ ∈ R, we obtain

c k = cϒ 2π k 2M -ϕ -ω 1 -cϒ 2π k 2M -ϕ -ω 2 . (52)
Moreover, we have ϕ = arg(e -jM ω 1 + e -j(M ω 2 +π ) )/M, so ϕ = -(ω 1 + ω 2 )/2 + 2π /(4M) + 2πk 0 /(2M), for some k 0 ∈ Z.

An example is illustrated in Figure 2. In Figure 3, we show the behavior of the atomic norm of , as a function of d(z 1 , z 2 ). We observe, without proof, that if this distance is larger than or equal to 1/(2M), the measure δ z 1δ z 2 is identifiable.

Optimal transport controlled in Fourier domain

Building upon the analysis in the previous section, we now focus on optimal transport between two measures, given only finite sequences of their Fourier coefficients. Optimal transport of k =1 c k a(z k ), so a = µ ⋆ TV and µ ⋆ is the measure of minimal total variation norm explaining .

measures [START_REF] Villani | Topics in Optimal Transportation[END_REF][START_REF] Santambrogio | Optimal Transport for Applied Mathematicians[END_REF] has grown significantly in recent years, due to a wide range of applications, from image processing [START_REF] Rabin | Wasserstein barycenter and its application to texture mixing[END_REF] to speech processing [START_REF] Alaux | Unsupervised hyper-alignment for multilingual word embeddings[END_REF].

Given a nonnegative cost function f on T 2 , the optimal transport cost T f (µ, µ ′ ) between two positive measures µ and µ ′ of M 2 , with same mass µ(T) = µ ′ (T), is the infimum, over all positive measures ν on T 2 with marginals µ and µ ′ on the first and second factors, respectively, of the cost ∫ T 2 f (z, z ′ )dν(z, z ′ ) [START_REF] Villani | Topics in Optimal Transportation[END_REF][START_REF] Rabin | Transportation distances on the circle[END_REF]. Let us extend this definition to signed measures. We define the optimal transport cost T f between a pair of signed measures (µ, µ ′ ) ∈ M 2 with same mass as the infimum, over all finite signed Borel measures ν on T 2 with marginals µ and µ ′ , of ∫ T 2 f (z, z ′ )d|ν |(z, z ′ ), where |ν | = ν + + ν -. This definition, however, is meaningful only for concave costs, since it vanishes for convex costs: if the cost is concave, i.e. f = ϕ • d for some increasing and concave function ϕ with ϕ(0) = 0, then T f (µ, µ ′ ) = T f (µ + + µ ′-, µ ′+ + µ -).

Some classical cost functions f are the following. First, if we consider the binary distance

f (z, z ′ ) = {0 if z = z ′ , 1 else}, then T f (µ, µ ′ ) is equal to the scaled Radon distance d R (µ, µ ′ ) defined as d R (µ, µ ′ ) = 1 2 µ -µ ′ TV , (53)
for every µ and µ ′ in M with same mass. Second, if f = d, for positive measures µ and µ ′ with same mass, then T f becomes the Wasserstein-1 distance, denoted by d W1 . Since f is a concave function of d, T f can be extended to signed measures, not only positive ones. Finally, if f = d 2 , for positive measures with same mass, the square root of T f is the Wasserstein-2 distance denoted by d W2 . We are not able to extend this functional to signed measures, however. k =1 c k δ z k given by Proposition 4, have the same sequence of Fourier coefficients = ( m ) M m=-M . In black, the atomic norm of . We observe that when d(z 1 , z 2 ) < 1/(2M), µ ⋆ is the measure of minimal total variation norm explaining , so a = µ ⋆ TV = 2 sin π Md(z 1 , z 2 ) . Else, when d(z 1 , z 2 ) ≥ 1/(2M), µ ♯ is the measure of minimal total variation norm explaining , so a = µ ♯ TV = 2.

Since we want to do optimal transport between measures, given only finite sequences of their Fourier coefficients, we now introduce the atomic transport cost, denoted by T a,f , associated to any nonnegative cost function f on T 2 . We restrict ourselves to the concave case; that is, f = ϕ • d for some increasing and concave function ϕ with ϕ(0) = 0. So, for every ( ,

′ ) ∈ V 2 with 0 = ′ 0 , T a,f ( , ′ ) := min T f (µ, µ ′ ) : (µ, µ ′ ) ∈ M 2 , Fµ = , Fµ ′ = ′ . ( 54
)
The difficulty is of course to exhibit computable expressions of the atomic transport cost.

In the simple case of the binary distance f (z, z ′ ) = {0 if z = z ′ , 1 else}, we call the atomic transport cost the atomic Radon distance, denoted d a0 , and we have, for every ( ,

′ ) ∈ V 2 with 0 = ′ 0 , d a0 ( , ′ ) := min d R (µ, µ ′ ) : (µ, µ ′ ) ∈ M 2 , Fµ = , Fµ ′ = ′ (55) = min 1 2 µ -µ ′ TV : (µ, µ ′ ) ∈ M 2 , Fµ = , Fµ ′ = ′ (56) = min 1 2 η TV : η ∈ M, Fη = -′ (57) = 1 2 -′ a . (58)
This distance, which is convex with respect to the pair ( , ′ ), allows us to quantify whether and ′ are different. Its precision naturally depends on M. We show in Figure 4 the value of d a0 taken at a pair of atoms, as a function of the distance of their frequencies, which is explained by the analysis carried out in Section 9. Let us now turn our attention to the 1-Wasserstein distance. When f = d, we call the atomic transport cost the atomic 1-Wasserstein distance, denoted by d a1 , defined as follows: for every

( , ′ ) ∈ V 2 with 0 = ′ 0 , d a1 ( , ′ ) := min d W1 (µ, µ ′ ) : (µ, µ ′ ) ∈ M 2 , Fµ = , Fµ ′ = ′ . ( 59
)
Like for the atomic Radon distance, we can express the atomic 1-Wasserstein distance as a semidefinite program:

Proposition 5. Let ( , ′ ) ∈ V 2 with 0 = ′ 0 . Set w = ( m -′ m )/(j2πm)
M m=-M , with w 0 = 0, and set W = T(w). Then

d a1 ( , ′ ) = min η TV : η ∈ M, Fη = w with (60) j2πmw m = m -′ m , m = -M, . . . , M , = min X ,β ∈R 2 M+1 tr(X ) + β s.t. X ∈ T(V) (61) and X 0 and X -W + βId 0, = min X 2 M+1 tr(X ) + i + (W -X ) s.t. (62)
X ∈ T(V) and X 0, where i + denotes the largest eigenvalue.

Proof: For every (µ, µ ′ ) ∈ M 2 , it is known that d W1 (µ, µ ′ ) only depends on µµ ′ and has the expression where F and F ′ are the cumulative distribution functions of µ and µ ′ , respectively [START_REF] Rabin | Transportation distances on the circle[END_REF][START_REF] Cabrelli | The Kantorovich metric for probability measures on the circle[END_REF]. Let α ∈ R and set the function η = F -F ′α. Set = Fµ and ′ = Fµ ′ . Let us interpret η as a measure in M and set w = Fη. We have the relationship j2πmw m = m -′ m , for every m = -M, . . . , M. The expression (60) follows. Then we see that d a1 ( , ′ ) takes a form similar to the atomic norm of w, but with w 0 , or equivalently the real α in (63) or the real β in (61), as an additional degree of freedom. So, ( 61) is nothing but the semidefinite characterization [START_REF] Auton | Investigation of procedures for automatic resonance extraction from noisy transient electromagnetics data[END_REF] of the atomic norm. Finally, for a fixed X , the real β minimizing β subject to X -W + βId 0 is the largest eigenvalue of W -X , so that we can remove this variable.

d W1 = min α ∈R 1 2π ∫ T |F (z) -F ′ (z) -α |dz, (
In Figure 5, we can see that when and ′ are atoms, d a1 ( , ′ ) closely matches the distance between their frequencies.

Finally, let us focus on the case of the squared 2-Wasserstein distance d 2 W2 = T f with f = d 2 . Since f is not concave with respect to d, the optimal transport cost T f (µ, µ ′ ) does not only depend on µµ ′ [START_REF] Rabin | Transportation distances on the circle[END_REF][START_REF] Cabrelli | The Kantorovich metric for probability measures on the circle[END_REF]. Moreover, it seems difficult to extend it to a convex functional with respect to a pair of signed measures. So, we restrict ourselves to the positive case. That is, given We can note that this distance is exact for atoms: for every z ∈ T, there is only one positive measure explaining a(z), which is δ z , and there is only one positive measure ν on T 2 having two Diracs as marginals, which is a Dirac. So, for every (z, z ′ ) ∈ T 2 , d a2 a(z), a(z ′ )

( , ′ ) ∈ V 2 with 0 = ′ 0 , T
2 = d(z, z ′ ) 2 .
However, there is little hope to express the squared atomic 2-Wasserstein distance as a finitedimensional convex program. So, we will propose an approximation of this functional. We restrict ourselves to the case where ′ is fixed and is an atom: ′ = a(z ′ ) for some z ′ ∈ T. We want to design a function da2,z ′ of ∈ V, with 0 = 1 and T( ) 0, whose square approximates d a2 ( , ′ ) 2 = min{ ∫ T d(z, z ′ ) 2 dµ(z) : µ is a positive measure and Fµ = }. For this, we mimic the construction of d a1 , which is the total variation of a measure whose 'first derivative' explains -′ , up to a constant: we define da2,z ′ ( ) 2 as the total variation of a measure whose 'second derivative' explains m - M m=-M , with w 0 = 0, and W = T(w); then da2,z ′ ( ) 2 = i + (-W ).

(66)

In Figure 6, we can see that when = a(z) is an atom, da2,z ′ ( ) 2 closely matches the squared distance d(z, z ′ ) 2 between the frequencies.

In Figure 7, we show a toy application of the proposed approximate atomic 2-Wasserstein distance: we compute the Wasserstein barycenter [START_REF] Agueh | Barycenters in the Wasserstein space[END_REF] ∈ V of two atoms a(z) and a(z ′ ), with z = e j2π /10 , z ′ = e j2π /20 , M = 10. That is, minimizes da2,z ( ) 2 + da2,z ′ ( ) 2 . The overrelaxed Chambolle-Pock algorithm is used again to solve this convex optimization problem. The solution turns out to be the atom a(e j3π /20 ), whose frequency is the mean of the frequencies of the two given atoms. So, we can interpolate between two sampled complex exponentials, by solving a convex optimization problem, without having to first identify their frequencies. This opens the door to many signal processing applications. parts of the values are in ordinate, with black and red colors, respectively. In (c): the atomic Wasserstein barycenter ∈ V, which minimizes da2,z ( ) 2 + da2,z ′ ( ) 2 ; it turns out that is the atom a(e j3π /20 ), whose frequency is the mean of the frequencies of the two given atoms.

Conclusion

We studied the properties of decompositions of vectors into complex exponentials and we proposed convex formulations to implement these decompositions. Roughly speaking, the analysis is based on splitting the vectors (and associated Toeplitz matrices) into negative and positive (semidefinite) parts and applying the powerful characterization of Toeplitz matrices by Carathéodory to them. Besides, we have introduced functionals to perform optimal transport between two measures, given only a finite number of their Fourier coefficients. This makes it possible to deal with measures computationally, without the need to discretize them on a grid of predefined locations. Further theoretical investigations are needed and the relationship between the proposed construction and existing spectral transport approaches for audio signal processing [START_REF] Flamary | Optimal spectral transportation with application to music transcription[END_REF][START_REF] Rolet | Blind source separation with optimal transport non-negative matrix factorization[END_REF] should be investigated. In future work, the author will apply this framework to the approximation and reconstruction of signals and images [START_REF] Condat | A convex approach to K-means clustering and image segmentation[END_REF][START_REF] Polisano | A convex approach to super-resolution and regularization of lines in images[END_REF][START_REF] Condat | A convex lifting approach to image phase unwrapping[END_REF].

Remark 9 .

 9 Let n + (T( )) and n -(T( )) be the number of positive and negative eigenvalues of T( ), respectively. Then every atomic decomposition = K k =1 c k a(z k ) has at least n + (T( )) positive c k and at least n -(T( )) negative c k .

µ

  TV s.t. m = ∫ T z -m dµ(z), m = -M, . . . , M, (31) has a solution, concentrated at a finite number K ≤ 2M of points of T, i.e. of the form µ ⋆ = K k =1 c k δ z k , for distinct z k ∈ T and nonzero reals c k . Equivalently, a minimal atomic decomposition = K k =1 c k a(z k ) exists. We obtain the parameters z k and c k by combining those of + and -

Figure 2 :

 2 Figure 2: The two measures µ ♯ = δ e j2π 0.51δ e j2π 0.54 , in red, and µ ⋆ = 2M k =1 c k δ z k given by Propositions 3 and 4, in blue, for M = 10, have the same sequence of Fourier coefficients = ( m ) M m=-M . These two measures have total variation norm µ ♯ TV = 2 and µ ⋆ TV = 2M k =1 |c k | ≈ 1.6, respectively. Numerical minimization of the atomic norm of yields the minimal atomic decomposition = 2Mk =1 c k a(z k ), so a = µ ⋆ TV and µ ⋆ is the measure

Figure 3 :

 3 Figure 3: Value of the total variation norm µ ♯ TV = {0 if z 1 = z 2 , 2 else}, in red, and µ ⋆ TV = 2M k =1 |c k |, in blue, as a function of d(z 1 , z 2 ), for M = 10, where the two measures µ ♯ = δ z 1δ z 2 and µ ⋆ = 2Mk =1 c k δ z k given by Proposition 4, have the same sequence of Fourier coefficients = ( m ) M m=-M . In black, the atomic norm of . We observe that when d(z 1 , z 2 ) < 1/(2M), µ ⋆ is the

  ( ) 0 and T( ) 0, we define the atomic 2-Wasserstein distance as d a2 ( , ′ ) := min d W2 (µ, µ ′ ) : µ and µ ′ are positive measures, (64) Fµ = , Fµ ′ = ′ .

Figure 6 :

 6 Figure 6: Value of the approximate squared atomic 2-Wasserstein distance da2,z 2 a(z 1 )2 as a function of d(z 1 , z 2 ), in black, with M = 10. We observe that the this function is very close to d(z 1 , z 2 ) 2 , in red, except for the largest values, where the squared distance is underestimated.

  ; that is, da2,z ′ ( ) 2 := min η(T) : η is a positive measure, Fη = w,with -4π 2 m 2 w m = m -2 ′ m + ( ′ m ) 2 * m , m = -M, . . . , M .(65)It turns out that this function has a closed-form expression: letw = ( m -2 ′ m + ( ′ m ) 2 * m )/(-4π 2 m 2 )

Figure 7 :

 7 Figure 7: In (a) and (b): two atoms a(z) = a(z) m M m=-M and a(z ′ ) = a(z ′ ) m M m=-M , respectively, with z = e j2π /10 , z ′ = e j2π /20 , M = 10; the index m is in abscissa, the real and imaginary

  If T( ) is indefinite, we have seen above that well separated measures are identifiable, but this is only a sufficient condition: there are measures which are not well separated but are identifiable. The picture is different if we restrict ourselves to positive measures; as a consequence of Proposition 2, we have the property: a positive measure is identifiable if and only if it is concentrated at K ≤ M points of T. It is remarkable that all the information about a positive discrete measure concentrated in at most M points is encoded in its coefficients ( m ) M m=1 (since 0 can be recovered as the opposite of the smallest eigenvalue of T((

* M , . . . , * 1 , 0, 1 , . . . , M )) ), characterized by the critical number 2M of real degrees of freedom, without any constraint.

  Figure4: Value of the atomic Radon distance d a0 a(z 1 ), a(z 2 ) as a function of d(z 1 , z 2 ), in black, with M = 10. We observe that if d(z 1 , z 2 ) < 1/(2M), d a0 a(z 1 ), a(z 2 ) = sin π Md(z 1 , z 2 ) . However, if d(z 1 , z 2 ) ≥ 1/(2M), d a0 a(z 1 ), a(z 2 ) = 1, which is the expected value, namely the Radon distance between δ z 1 and δ z 2 .

	1
	0.8
	0.6
	0.4
	0.2
	0
	0	0.05	0.1	0.15	0.2	0.25	0.3	0.35	0.4	0.45	0.5
						d(z 1 ,z 2 )

  Figure5: Value of the atomic 1-Wasserstein distance d a1 a(z 1 ), a(z 2 ) as a function of d(z 1 , z 2 ), in black, with M = 10. We observe that the atomic 1-Wasserstein distance is equal to d(z 1 , z 2 ), in red, except for the largest values, where this distance is underestimated.

	0.5
	0.45
	0.4
	0.35
	0.3
	0.25
	0.2
	0.15
	0.1
	0.05
	0
	0	0.05	0.1	0.15	0.2	0.25	0.3	0.35	0.4	0.45	0.5
						d(z 1 ,z 2 )
											63)