Claudio Durastanti

Domenico Marinucci
email: marinucc@mat.uniroma2.it.

Giovanni Peccati

Normal Approximations for Wavelet Coefficients on Spherical Poisson Fields

Keywords: Berry-Esseen Bounds, Malliavin Calculus, Multidimensional Normal Approximation, Poisson Process, Stein's Method, Spherical Wavelets AMS Classification: 60F05, 42C40, 33C55, 60G60, 62E20

We compute explicit upper bounds on the distance between the law of a multivariate Gaussian distribution and the joint law of wavelets/needlets coefficients based on a homogeneous spherical Poisson field. In particular, we develop some results from Peccati and Zheng (2011), based on Malliavin calculus and Stein's methods, to assess the rate of convergence to Gaussianity for a triangular array of needlet coefficients with growing dimensions. Our results are motivated by astrophysical and cosmological applications, in particular related to the search for point sources in Cosmic Rays data.

Introduction

The aim of this paper is to establish multidimensional normal approximation results for vectors of random variables having the form of wavelet coefficients integrated with respect to a Poisson measure on the unit sphere. The specificity of our analysis is that we require the dimension of such vectors to grow to infinity. Our techniques are based on recently obtained bounds for the normal approximation of functionals of general Poisson measures (see [START_REF] Peccati | Stein's method and normal approximation of Poisson functionals[END_REF][START_REF] Peccati | Multi-dimensional Gaussian fluctuations on the Poisson space[END_REF]), as well as on the use of the localization properties of wavelets systems on the sphere (see [START_REF] Narcowich | Localized Tight Frames on Spheres[END_REF], as well as the recent monograph [START_REF] Marinucci | Random Fields on the Sphere. Representation, Limit Theorems and Cosmological Applications[END_REF]). A large part of the paper is devoted to the explicit determination of the above quoted bounds in terms of dimension.

Motivation and overview

A classical problem in asymptotic statistics is the assessment of the speed of convergence to Gaussianity (that is, the computation of explicit Berry-Esseen bounds) for parametric and nonparametric estimation procedures -for recent references connected to the main topic of the present paper, see for instance [START_REF] Huang | Berry-Esseen bounds for kernel estimates of stationary processes[END_REF][START_REF] Li | Berry-Esseen bounds for wavelet estimator in a regression model with linear process errors[END_REF][START_REF] Yang | The Berry-Esséen type bound of sample quantiles for strong mixing sequence[END_REF]. In this area, an important novel development is given by the derivation of effective Berry-Esseen bounds by means of the combination of two probabilistic techniques, namely the Malliavin calculus of variations and the Stein's method for probabilistic approximations. The monograph [START_REF] Chen | Normal Approximation by Stein's Method[END_REF] is the standard modern reference for Stein's method, whereas [START_REF] Nourdin | Normal approximations using Malliavin calculus: from Stein's method to universality[END_REF] provides an exhaustive discussion of the use of Malliavin calculus for proving normal approximation results on a Gaussian space. The fact that one can use Malliavin calculus to deduce normal approximation bounds (in total variation) for functionals of Gaussian fields was first exploited in [START_REF] Nourdin | Stein's method on Wiener chaos[END_REF] -where one can find several quantitative versions of the "fourth moment theorem" for chaotic random variables proved in [START_REF] Nualart | Central limit theorems for sequences of multiple stochastic integrals[END_REF]. Lower bounds can also be computed, entailing that the rates of convergence provided by these techniques are sharp in many instances -see again [START_REF] Nourdin | Normal approximations using Malliavin calculus: from Stein's method to universality[END_REF].

In a recent series of contributions, the interaction between Stein's method and Malliavin calculus has been further exploited for dealing with the normal approximation of functionals of a general Poisson random measure. The most general abstract results appear in [START_REF] Peccati | Stein's method and normal approximation of Poisson functionals[END_REF] (for one-dimensional normal approximations) and [START_REF] Peccati | Multi-dimensional Gaussian fluctuations on the Poisson space[END_REF] (for normal approximations in arbitrary dimensions). These findings have recently found a wide range of applications in the field of stochastic geometry -see [START_REF] Lachièze-Rey | Fine Gaussian fluctuations on the Poisson space I: contractions, cumulants and geometric random graphs[END_REF][START_REF] Lachièze-Rey | Fine Gaussian fluctuations on the Poisson space II: rescaled kernels, marked processes and geometric U-statistics[END_REF][START_REF] Minh | Malliavin-Stein method for multi-dimensional U-statistics of Poisson point processes[END_REF][START_REF] Last | Moments and central limit theorems for some multivariate Poisson functionals[END_REF][START_REF] Reitzner | Central Limit Theorems for U-Statistics of Poisson Point Processes[END_REF] for a sample of geometric applications, as well as the webpage http://www.iecn.u-nancy.fr/∼nourdin/steinmalliavin.htm for a constantly updated resource on the subject.

The purpose of this paper is to apply and extend the main findings of [START_REF] Peccati | Stein's method and normal approximation of Poisson functionals[END_REF][START_REF] Peccati | Multi-dimensional Gaussian fluctuations on the Poisson space[END_REF] in order to study the multidimensional normal approximation of the elements of the first Wiener chaos of a given Poisson measure. Our main goal is to deduce bounds that are well-adapted to deal with applications where the dimension of a given statistic increases with the number of observations. This is a framework which arises naturally in many relevant fields of modern statistical analysis; in particular, our principal motivation originates from the implementation of wavelet systems on the sphere. In these circumstances, when more and more data become available, a higher number of wavelet coefficients is evaluated, as it is customarily the case when considering, for instance, thresholding nonparametric estimators. We shall hence be concerned with sequences of Poisson fields, whose intensity grows monotonically. We then exploit the wavelets localization properties to establish bounds that grow linearly with the number of functionals considered; we are then able to provide explicit recipes, for instance, for the number of joint testing procedures that can be simultaneously entertained ensuring that the Gaussian approximation may still be shown to hold, in a suitable sense.

Main contributions

Consider a sequence {X i : i ≥ 1} with values in the unit sphere S 2 , and define {ψ jk } to be the collection of the spherical needlets associated with a certain constant B > 1, see Section 3.1 below for more details and discussion. Write also σ 2 jk = E[ψ jk (X 1 ) 2 ] and b jk = E[ψ jk (X 1 )], and consider an independent (possibly inhomogeneous) Poisson process {N t : t ≥ 0} on the real line such that E[N t ] = R(t) → ∞, as t → ∞. Formally, our principal aim is to establish conditions on the sequences {j(n

) : n ≥ 1}, {R(n) : n ≥ 1} and {d(n) : n ≥ 1} ensuring that the distribution of the centered d(n)-dimensional vector Y n = (Y n,1 , ..., Y n,d(n) ) (1.1) = 1 R(n)   N (n) i=1 ψ j(n)k 1 (X i ) σ j(n)k 1 - R(n)b j(n)k 1 σ j(n)k 1 , ..., N (n) i=1 ψ j(n)k d(n) (X i ) σ j(n)k 1 - R(n)b j(n)k d(n) σ j(n)k d(n)  
is asymptotically close, in the sense of some smooth distance denoted d 2 (see Definition 2.6), to the law of a d(n)-dimensional Gaussian vector, say Z n , with centered and independent components having unit variance. The use of a smooth distance allows one to deduce minimal conditions for this kind of asymptotic Gaussianity. The crucial point is that we allow the dimension d(n) to grow to infinity, so that our results require to explicitly assess the dependence of each bound on the dimension. We shall perform our tasks through the following main steps: (i) Proposition 4.1 deals with one-dimensional normal approximations, (ii) Proposition 5.4 deals with normal approximations in a fixed dimension, and finally (iii) in Theorem 5.5 we deduce a bound that is well-adapted to the case d(n) → ∞. More precisely, Theorem 5.5 contains an upper bound linear in d(n), that is, an estimate of the type

d 2 (Y n , Z n ) ≤ C(n) × d(n) . (1.
2)

It will be shown in Corollary 5.6, that the sequence C(n) can be chosen to be n) can be viewed as a measure of the "effective sample size" for the components of Y n .

O 1/ R(n)B -2j(n) ; as discussed below in Remark 4.3, R(n) × B -2j(

About de-Poissonization

Our results can be used in order to deduce the asymptotic normality of de-Poissonized linear statistics with growing dimension. To illustrate this point, assume that the random variables X i are uniformly distributed on the sphere. Then, it is well known that b jk = 0, whenever j > 1. In this framework, when j(n) > 1 for every n, R(n) = n and d(n)/n 1/4 → 0, the conditions implying that Y n is asymptotically close to Gaussian, automatically ensure that the law of the de-Poissonized vector

Y ′ n = (Y ′ n,1 , ..., Y ′ n,d(n) ) = 1 √ n n k=1 ψ j(n)k 1 (X i ) σ j(n)k 1 , ..., n k=1 ψ j(n)k d(n) (X i ) σ j(n)k d(n) (1.3)
is also asymptotically close to Gaussian. The reason for this phenomenon is nested in the statement of the forthcoming (elementary) Lemma 1.1.

Lemma 1.1 Assume that R(n) = n, that the X i 's are uniformly distributed on the sphere, and that j(n) > 1 for every n. Then, there exists a universal constant M such that, for every n and every Lipschitz function ϕ : R d(n) → R, the following estimate holds:

E[ϕ(Y ′ n )] -E[ϕ(Y n )] ≤ M ϕ Lip d(n) n 1/4 .
Proof. Fix l = 1, ..., d(n), and write

β l (x) = ψ j(n)k l (x) σ j(n)k l , in such a way that E[β l (X 1 ) 2 ] = 1. One has that E[(Y ′ n,l -Y n,l ) 2 ] = 2(1 -α n ),
where

α n = 1 n n m=0 e -n n m m! (n ∧ m) = 1 - e -n n n n! .
This gives the estimate

E[|Y ′ n,l -Y n,l |] ≤ E[(Y ′ n,l -Y n,l ) 2 ] ≤ 2 e -n n n n! ,
so that the conclusion follows from an application of Stirling's formula and of the Lipschitz property of ϕ.

Remark 1.2 (i) Lemma 1.1 implies that one can obtain an inequality similar to (1.2) for Y ′ n , that is:

d 2 (Y ′ n , Z n ) ≤ C(n) + M n 1/4 × d(n).
(ii) With some extra work, one can obtain estimates similar to those in Lemma 1.1 also when the constants b j(n)k l are possibly different from zero. This point, that requires some lenghty technical considerations, falls slightly outside the scope of this paper and will be pursued in full generality elsewhere.

(iii) In [START_REF] Bentkus | On the dependence of the Berry-Esseen bound on dimension[END_REF], Bentkus proved the following (yet unsurpassed) bound. Assume that {X i : i ≥ 1} is a collection of i.i.d. d-dimensional vectors, such that X 1 is centered and with covariance equal to the identity matrix. Set

S n = n -1/2 (X 1 + • • • X n ),
C ⊂ R d E[1 C (S n )] -E[1 C (Z)] ≤ d 1/4 400β √ n , where β = E[ X 1 3 R d ].
It is unclear whether one can effectively use this bound in order to investigate the asymptotic Gaussianity of sequences of random vectors of the type (1.1)-(1.3), in particular because, for a fixed n, the components of Y n , Y ′ n have in general a non trivial correlation. Note also that a simple application of Jensen inequality shows that βd 1/4 n -1/2 ≥ d 7/4 n -1/2 . However, a direct comparison of Bentkus' estimates with our "linear" rate in d (see (1.2), as well as Theorem 5.5 below) is unfeasible, due to the differences with our setting, namely concerning the choice of distance, the structure of the considered covariance matrices, the Poissonized environment, and the role of B j(n) discussed in Remark 4.3 .

(iv) A careful inspection of the proofs of our main results reveals that the findings of this paper have a much more general validity, and in particular can be extended to kernel estimators on compact spaces satisfying mild concentration and equispacing properties (see also [START_REF] Jupp | Sobolev tests of goodness of fit of distributions on compact Riemannian manifolds[END_REF][START_REF] Jupp | Data-driven Sobolev tests of uniformity on compact Riemannian manifolds[END_REF]). In this paper, however, we decided to stick to the presentation on the sphere for definiteness, and to make the connection with applications clearer. Some more general frameworks are discussed briefly at the end of Section 5.

(v) For notational simplicity, throughout this paper we will stick to the case where all the components in our vector statistics are evaluated at the same scale j(n) (see below for more precise definitions and detailed discussion). The relaxation of this assumption to cover multiple scales (j 1 (n), ...j d (n)) does not require any new ideas and is not considered here for brevity's sake.

Plan

The plan of the paper is as follows: in Section 2 we provide some background material on Stein-Malliavin bounds in the case of Poisson random fields, and we describe a suitable setting for the current paper, entailing sequences of fields with monotonically increasing governing measures. We provide also some new results, ensuring that the Central Limit Theorems we are going to establish are stable, in the classical sense. In Section 3 we recall some background material on the construction of tight wavelet systems on the sphere (see [START_REF] Narcowich | Localized Tight Frames on Spheres[END_REF][START_REF] Narcowich | Decomposition of Besov and Triebel-Lizorkin Spaces on the Sphere[END_REF] for the original references, as well as [START_REF] Marinucci | Random Fields on the Sphere. Representation, Limit Theorems and Cosmological Applications[END_REF]Chapter 10]) and we explain how to express the corresponding wavelet coefficients in terms of stochastic integrals with respect to a Poisson random measure. We also illustrate shortly some possible statistical applications.

In Section 4 we provide our bounds in the one-dimensional case; these are simple results which could have been established by many alternative techniques, but still they provide some interesting insights into the "effective area of influence" of a single component of the wavelet system. The core of the paper is in Section 5, where the bound is provided in the multidimensional case, allowing in particular for the number of coefficients to be evaluated to grow with the number of observations. This result requires a careful evaluation of the upper bound, which is made possible by the localization properties in real space of the wavelet construction.

Poisson Random Measures and Stein-Malliavin Bounds

In order to study the asymptotic behaviour of linear functionals of Poisson measures on the sphere S 2 , we start by recalling the definition of a Poisson random measure -for more details, see for instance [START_REF] Peccati | Wiener chaos: moments, cumulants and diagrams[END_REF][START_REF] Privault | Stochastic analysis in discrete and continuous settings with normal martingales[END_REF][START_REF] Schneider | Stochastic and integral geometry[END_REF]. We work on a probability space (Ω, F , P ).

Definition 2.1 Let (Θ, A, µ) be a σ-finite measure space, and assume that µ has no atoms (that is, µ({x}) = 0, for every x ∈ Θ). A collection of random variables {N (A) : A ∈ A} , taking values in Z + ∪ {+∞} , is called a Poisson random measure (PRM) on Θ with intensity measure (or control measure) µ if the following two properties hold:

1. For every A ∈ A, N (A) has Poisson distribution with mean µ (A);

2. If A 1 , . . . A n ∈ A are pairwise disjoint, then N (A 1 ) , . . . , N (A n ) are independent.
Remark 2.2 (i) In Definition 2.1, a Poisson random variable with parameter λ = ∞ is implicitly set to be equal to ∞.

(ii) Points 1 and 2 in Definition 2.1 imply that, for every ω ∈ Ω, the mapping A → N (A, ω) is a measure on Θ. Moreover, since µ is non atomic, one has that

P N({x}) = 0 or 1, ∀x ∈ Θ = 1. (2.4)
Assumption 2.3 Our framework for the rest of the paper will be the following special case of Definition 2.1:

(a) We take Θ = R + × S 2 , with A = B(Θ), the class of Borel subsets of Θ.

(b) The symbol N indicates a Poisson random measure on Θ, with homogeneous intensity given by µ = ρ × ν, where ρ is some measure on R + and ν is a probability on S 2 of the form ν(dx) = f (x)dx, where f is a density on the sphere. We shall assume that ρ({0}) = 0 and that the mapping ρ → ρ([0, t]) is strictly increasing and diverging to infinity as t → ∞. We also adopt the notation

R t := ρ([0, t]), t ≥ 0, (2.5) that is, t → R t is the distribution function of ρ.
Remark 2.4

(i) For a fixed t > 0, the mapping

A → N t (A) := N([0, t] × A) (2.6)
defines a Poisson random measure on S 2 , with non-atomic intensity

µ t (dx) = R t • ν(dx) = R t • f (x)dx. (2.7)
Throughout this paper, we shall assume f (x) to be bounded and bounded away from zero, e.g.

ζ 1 ≤ f (x) ≤ ζ 2 , some ζ 1 , ζ 2 > 0 , for all x ∈ S 2 . (2.8)
(ii) Let {X i = i ≥ 1} be a sequence of i.i.d. random variables with values in S 2 and common distribution equal to ν. Then, for a fixed t > 0, the random measure

A → N t (A) = N([0, t] × A) has the same distribution as A → N i=1 δ X i (A)
, were δ x indicates a Dirac mass at x, and N is an independent Poisson random variable with parameter R t . This holds because: (a) since ν is a probability measure, the support of the random measure N t (written supp(N t )) is almost surely a finite set, and (b) conditionally on the event {N t (S 2 ) = n} (which is the same as the event { |supp(N t )| = n} -recall (2.4)), the points in the support of N t are distributed as n i.i.d. random variables with common distribution ν.

(iii) By definition, for every t 1 < t 2 one has that a random variable of the type

N t 2 (A) - N t 1 (A), A ⊂ S 2
, is independent of the random measure N t 1 , as defined in (2.6).

(iv) To simplify the discussion, one can assume that ρ(ds) = R • ℓ(ds), where ℓ is the Lebesgue measure and R > 0, in such a way that R t = R • t.

We will now introduce two distances between laws of random variables taking values in R d . Both distances define topologies, over the class of probability distributions on R d , that are strictly stronger than convergence in law. One should observe that, in this paper, the first one (Wasserstein distance) will be only used for random elements with values in R. Given a function g ∈ C 1 (R d ), we write g Lip = sup

x∈R d ∇g(x) R d . If g ∈ C 2 (R d ), we set M 2 (g) = sup
X, Y with values in R d (d ≥ 1) and such that E X R d , E Y R d < ∞, is given by: d W (X, Y ) = sup g: g Lip ≤1 |E [g (X)] -E [g (Y )]| , Definition 2.6 The distance d 2 between the laws of two random vectors X, Y with values in R d (d ≥ 1), such that E X R d , E Y R d < ∞, is given by: d 2 (X, Y ) = sup g∈H |E [g (X)] -E [g (Y )]| , where H denotes the collection of all functions g ∈ C 2 R d such that g Lip ≤ 1 and M 2 (g) ≤ 1.
We now present, in a form adapted to our goals, two upper bounds involving random variables living in the so-called first Wiener chaos of N. The first bound was proved in [START_REF] Peccati | Stein's method and normal approximation of Poisson functionals[END_REF], and concerns normal approximations in dimension 1 with respect to the Wasserstein distance. The second bound appears in [START_REF] Peccati | Multi-dimensional Gaussian fluctuations on the Poisson space[END_REF], and provides estimates for multidimensional normal approximations with respect to the distance d 2 . Both bounds are obtained by means of a combination of the Malliavin calculus of variations and the Stein's method for probabilistic approximations.

Remark 2.7 (i) Let f ∈ L 2 (Θ, µ) ∩ L 1 (Θ, µ).
In what follows, we shall use the symbols N(f ) and N (f ), respectively, to denote the Wiener-Itô integrals of f with respect to N and with respect to the compensated Poisson measure

N(A) = N(A) -µ(A), A ∈ B(Θ), (2.9)
where one uses the convention

N(A) -µ(A) = ∞ whenever µ(A) = ∞ (recall that µ is σ-finite).
Note that, for N(f ) to be well-defined, one needs that f ∈ L 1 (Θ, µ), whereas for N (f ) to be well-defined one needs that f ∈ L 2 (Θ, µ). We will also make use of the following isometric property: for every f, g ∈L 2 (Θ, µ),

E[ N(f ) N (g)] = Θ f (x)g(x)µ(dx). (2.10)
The reader is referred e.g. to [START_REF] Peccati | Wiener chaos: moments, cumulants and diagrams[END_REF]Chapter 5] for an introduction to Wiener-Itô integrals.

(ii) For most of this paper, we shall consider Wiener-Itô integrals of functions f having the form f = [0, t] × h, where t > 0 and h ∈ L 2 (S 2 , ν) ∩ L 1 (S 2 , ν). For a function f of this type one simply writes

N(f ) = N([0, t] × h) := N t (h), and N(f ) = N ([0, t] × h) := Nt (h). (2.11)
Observe that this notation is consistent with the one introduced in (2.6). Indeed, it is easily seen that N t (h) (resp. Nt (h)) coincide with the Wiener-Itô integral of h with respect to N t (resp. with respect to the compensated measure

Nt = N t -µ t = N t -R t • ν). (iii) In view of Remark 2.4-(ii), one also has that, for h ∈ L 2 (S 2 , ν) ∩ L 1 (S 2 , ν), N t (h) = x∈supp(Nt) h(x), and
Nt (h) = x∈supp(Nt) h(x) - S 2 h(x)µ t (dx), (2.12)
with µ t defined as in (2.7).

Theorem 2.8 Let the notation and assumptions of this section prevail.

1. Let h ∈ L 2 (S 2 , ν) := L 2 (ν), let Z ∼ N (0, 1
) and fix t > 0. Then, the following bound holds (remember the definition (2.7)):

d W ( Nt (h), Z) ≤ 1 -h 2 L 2 (S 2 ,µt) + S 2 |h(z)| 3 µ t (dz). (2.13) As a consequence, if {h t } ⊂ L 2 (ν) ∩ L 3 (ν) is a collection of kernels verifying, as t → ∞, h t L 2 (S 2 ,µt) → 1 and h t L 3 (S 2 ,µt) → 0, (2.14)
one has the CLT

N(h t ) Law -→ Z, (2.15)
and the inequality (2.13) provides an explicit upper bound in the Wasserstein distance.

2. For a fixed integer d ≥ 1, let Y ∼ N d (0, C), with C positive definite and let

F t = (F t,1 , . . . , F t,d ) = Nt (h t,1 ) , . . . Nt (h t,d )
be a collection of d-dimensional random vectors such that h t,a ∈ L 2 (ν). If we call Γ t the covariance matrix of F t , that is,

Γ t (a, b) = E Nt (h t,a ) Nt (h t,b ) = h t,a , h t,b L 2 (S 2 ,µt) , a, b = 1, ..., d, then: d 2 (F t , Y ) ≤ C -1 op C 1 2 op C -Γ t H.S.
(2.16)

+ √ 2π 8 C -1 3 2 op C op d i,j,k=1 S 2 |h t,i (x)| |h t,j (x)| |h t,k (x)| µ t (dx) , ≤ C -1 op C 1 2 op C -Γ t H.S.
(2.17)

+ d 2 √ 2π 8 C -1 3 2 op C op d i=1 S 2 |h t,i (x)| 3 µ t (dx) ,
where • op and • H.S. stand, respectively, for the operator and Hilbert-Schmidt norms. In particular, if

Γ t (a, b) -→ C (a, b) and S 2 |h t,a (x)| 3 µ t (dx) -→ 0 as t -→ ∞, for a, b = 1, . . . d, then d 2 (F t , Y ) -→ 0 and F t converges in distribution to Y .
Remark 2.9 The estimate (2.16) will be used to deduce one of the main multidimensional bounds in the present paper. It is a direct consequence of Theorem 3.3 in [START_REF] Peccati | Multi-dimensional Gaussian fluctuations on the Poisson space[END_REF], where the following relation is proved: for every vector (F 1 , ..., F d ) of sufficiently regular centered functionals of Nt ,

d 2 (F, X) ≤ C -1 op C 1/2 op d i,j E C(i, j) -DF i , -DL -1 F j L 2 (µt) 2 + √ 2π 8 C -1 3/2 op C op S 2 µ t (dz)E   d i=1 |D z F i | 2 d j=1 D z L -1 F j   ,
where

D z F (ω) = F z (ω) -F (ω) , a.e. -µ(dz)P (dω) ,
and

F z (N) = F z (N + δ z ),
that is, the random variable F z is obtained by adding to the argument of F (which is a function of the point measure N), a Dirac mass at z, and L -1 is the so-called pseudo-inverse of the Ornstein-Uhlenbeck operator. The estimate (2.16) is then obtained by observing that, when

F i = F t,i = Nt (h t,i ), then D z F i = -D z L -1 F = h t,i (z), in such a way that d i,j E C(i, j) -DF i , -DL -1 F j L 2 (µ) 2 = C -K t H.S. , and
S 2 µ t (dz)E   d i=1 |D z F i | 2 d j=1 D z L -1 F j   = d i,j,k=1 S 2 |h t,i (x)| |h t,j (x)| |h t,k (x)| µ t (dx) .
The next statement deals with the interesting fact that the convergence in law implied by Theorem 2.8 is indeed stable, as defined e.g. in the classic reference [START_REF] Jacod | Limit Theorems for Stochastic Processes[END_REF]Chapter 4].

Proposition 2.10 The central limit theorem described at the end of Point 2 of Theorem 2.8 (and a fortiori the CLT at Point 1 of the same theorem) is stable with respect to σ(N) (the σ-field generated by N) in the following sense: for every random variable X that is σ(N)-measurable, one has that

(X, F t ) Law -→ (X, Y ), where Y ∼ N d (0, C) is independent of N.
Proof. We just deal with the case d = 1, the extension to a general d following from elementary considerations. A density argument shows that it is enough to prove the following claim:

if N(h n ) (h n ∈ L 2 (µ), n ≥ 1) is a sequence of random variables verifying E[ N (h n ) 2 ] = h n 2 L 2 (µ) → 1 and Θ |h n | 3 dµ → 0,
then for every fixed f ∈ L 2 (µ), the pair ( N(f ), N(h n )) converges in distribution, as n → ∞, to ( N(f ), Z), where Z ∼ N (0, 1) is independent of N. To see this, we start with the explicit formula (see e.g. [41, formula (5.3.31)]): for every λ, γ ∈ R

ψ n (λ, γ) := E[exp(iλ N (f ) + γ N(h n ))] = exp Θ e iλf (x)+iγhn(x) -1-i(λf (x) + γh n (x)) µ(dx) .
Our aim is to prove that, under the stated assumptions,

lim n→∞ log(ψ n (λ, γ)) = Θ e iλf (x) -1 -iλf (x) µ(dx) - γ 2 2 . Standard computations show that log(ψ n (λ, γ)) - Θ e iλf (x) -1 -iλf (x) µ(dx) - γ 2 2 ≤ γ 2 2 - γ 2 2 Θ h n (x) 2 µ(dx) + |γλ| | h n , f L 2 (µ) | + |γ| 3 6 Θ |h n (x)| 3 µ(dx) .
Since Θ |h n (x)| 3 µ(dx) → 0 and the mapping n → h n 2 L 2 (µ) is bounded, one has that h n , f L 2 (µ) → 0, and the conclusion follows by using the fact that h n 2 L 2 (µ) → 1 by assumption.

3 Needlet coefficients

Background: the needlet construction

We now provide an overview of the construction of the set of needlets on the unit sphere. The reader is referred to [START_REF] Marinucci | Random Fields on the Sphere. Representation, Limit Theorems and Cosmological Applications[END_REF]Chapter 10] for an introduction to this topic. Relevant references on this subject are: the seminal papers [START_REF] Narcowich | Localized Tight Frames on Spheres[END_REF][START_REF] Narcowich | Decomposition of Besov and Triebel-Lizorkin Spaces on the Sphere[END_REF], where needlets have been first defined; [START_REF] Geller | Continuous Wavelets on Manifolds[END_REF][START_REF] Geller | Nearly Tight Frames and Space-Frequency Analysis on Compact Manifolds[END_REF][START_REF] Geller | Spin Wavelets on the Sphere[END_REF][START_REF] Geller | Band-Limited Localized Parseval Frames and Besov Spaces on Compact Homogeneous Manifolds[END_REF], among others, for generalizations to homogeneous spaces of compact groups and spin fiber bundles; [START_REF] Baldi | Asymptotics for Spherical Needlets[END_REF][START_REF] Baldi | Subsampling Needlet Coefficients on the Sphere[END_REF][START_REF] Lan | On the Dependence Structure of Wavelet Coefficients for Spherical Random Fields[END_REF][START_REF] Mayeli | Asymptotic Uncorrelation for Mexican Needlets[END_REF] for the analysis of needlets on spherical Gaussian fields, and [START_REF] Marinucci | Spherical Needlets for CMB Data Analysis[END_REF][START_REF] Pietrobon | Integrated Sachs-Wolfe Effect from the Cross Correlation of WMAP3 Year and the NRAO VLA Sky Survey Data: New Results and Constraints on Dark Energy[END_REF][START_REF] Delabrouille | A Full Sky, Low Foreground, High Resolution CMB Map from WMAP[END_REF][START_REF] Feeney | First observational tests of eternal inflation: Analysis methods and WMAP 7-year results[END_REF] for some (among many) applications to cosmological and astrophysical issues; see also [START_REF] Mcewen | Cosmological applications of a wavelet analysis on the sphere[END_REF][START_REF] Roşca | Wavelet bases on the sphere obtained by radial projection[END_REF] for other approaches to spherical wavelets construction.

(Spherical harmonics) In Fourier analysis, the set of spherical harmonics {Y lm : l ≥ 0, m = -l, ..., l} provides an orthonormal basis for the space of square-integrable functions on the unit sphere L 2 (S 2 , dx) := L 2 (S 2 ), where dx stands for the Lebesgue measure on S 2 (see for instance [START_REF] Adler | Random Fields and Geometry[END_REF][START_REF] Koo | Sharp Adaptation for Spherical Inverse Problems with Applications to Medical Imaging[END_REF][START_REF] Marinucci | Random Fields on the Sphere. Representation, Limit Theorems and Cosmological Applications[END_REF][START_REF] Stein | Introduction to Fourier Analysis on Euclidean Spaces[END_REF]). Spherical harmonics are defined as the eigenfunctions of the spherical Laplacian ∆ S 2 corresponding to eigenvalues -l (l + 1), e.g. ∆ S 2 Y lm = -l(l + 1)Y lm , see again [START_REF] Marinucci | Random Fields on the Sphere. Representation, Limit Theorems and Cosmological Applications[END_REF][START_REF] Stein | Introduction to Fourier Analysis on Euclidean Spaces[END_REF][START_REF] Varshalovich | Quantum Theory of Angular Momentum[END_REF] for analytic expressions and more details and properties. For every l ≥ 0, we define K l as the linear space given by the restriction to the sphere of the polynomials with degree at most l. Plainly, one has that

K l = l k=0 span {Y km : m = -k, ..., k} ,
where the direct sum is in the sense of L 2 (S 2 ).

(Cubature points) It is well-known that for every integer l = 1, 2, ... there exists a finite set of cubature points Q l ⊂ S 2 , as well as a collection of weights {λ η }, indexed by the elements of Q l , such that

∀f ∈ K l , S 2 f (x)dx = η∈Q l λ η f (η).

Now fix B > 1, and write [x] to indicate the integer part of a given real x.

In what follows, we shall denote by X j = {ξ jk } and {λ jk }, respectively, the set Q [2B j+1 ] and the associated class of weights. We also write K j = card{X j }. As proved in [START_REF] Narcowich | Localized Tight Frames on Spheres[END_REF][START_REF] Narcowich | Decomposition of Besov and Triebel-Lizorkin Spaces on the Sphere[END_REF], cubature points and weights can be chosen to satisfy

λ jk ≈ B -2j , K j ≈ B 2j , ( 3.18)
where by a ≈ b, we mean that there exists c 1 , c 2 > 0 such that c 1 a ≤ b ≤ c 2 a (see also e.g. [START_REF] Baldi | Adaptive Density Estimation for directional Data Using Needlets[END_REF][START_REF] Pesenson | Sampling of band-limited vectors[END_REF][START_REF] Pesenson | Poincaré-type inequalities and reconstruction of Paley-Wiener functions on manifolds[END_REF] and [START_REF] Marinucci | Random Fields on the Sphere. Representation, Limit Theorems and Cosmological Applications[END_REF]Chapter 10]).

(Spherical needlets) Fix B > 1 as before, as well as a real-valued mapping b on (0, ∞).

We assume that b verifies the following properties: (i) the function b (•) has compact support in [B -1 , B] (in such a way that the mapping l → b l B j has compact support in l ∈ [B j-1 , B j+1 ]) (ii) for every ξ ≥ 1, ∞ j=0 b 2 (ξB -j ) = 1 (partition of unit property), and (iii) b (•) ∈ C ∞ (0, ∞). The collection of spherical needlets {ψ jk }, associated with B and b(•), are then defined as a weighted convolution of the projection operator L l ( x, y ) = l m=-l Y lm (x) Y lm (y), that is

ψ jk (x) := λ jk l b l B j L l ( x, ξ jk ) , (3.19)
(Localization) The properties of b entail the following quasi-exponential localization property (see [START_REF] Narcowich | Localized Tight Frames on Spheres[END_REF] or [30, Section 13.3]): for any τ = 1, 2, ... there exists κ τ > 0 such that for any x ∈ S 2 ,

|ψ jk (x)| ≤ κ τ B j (1 + B j arccos ( x, ξ jk )) τ , ( 3.20)
where d(x, y) := arccos ( x, y ) is the spherical distance. From localization, the following bound can be established on the L p (S 2 ) norms: for all 1 ≤ p ≤ +∞, there exist two positive constants q p and q ′ p such that

q p B j(1-2 p ) ≤ ψ jk Lp(S 2 ) ≤ q ′ p B j(1-2 p ) . (3.21)
(Needlets as frames) Finally, the fact that b is a partition of unit, allows on to deduce the following reconstruction formula (see again [START_REF] Narcowich | Localized Tight Frames on Spheres[END_REF]): for f ∈ L 2 (S 2 ):

f (x) = j,k β jk ψ jk (x) ,
where the convergence of the series is in L 2 (S 2 ), and

β jk := f, ψ jk L 2 (S 2 ) = S 2 f (x) ψ jk (x) dx , (3.22)
represents the so-called needlet coefficient of index j, k.

Two motivations: density estimates and point sources

The principal aim of this paper is to establish multidimensional asymptotic results for some possibly randomized version of random variables of the type

β jk = β (n) jk = 1 n n i=1 ψ jk (X i ) , j = 1, 2, ..., k = 1, ..., K j , (3.23)
where the function ψ jk is defined according to (3.19), and {X i : i ≥ 1} is some adequate sequence of i.i.d. random variables. We may also study the asymptotic behaviour, as t → ∞, of multi-dimensional object of the type β jk , k = 1, 2, ..., K j (t) , where t → K j (t) is a non-decreasing mapping possibly diverging to infinity, and j may change with t. In other words, as happens in realistic experimental circumstances, we may decide to focus on a growing number of coefficients as the number of (expected) events increase. Two strong motivations for this analysis, both coming from statistical applications, are detailed below.

(Density estimates) Consider a density function f on the sphere S 2 , that is: f is a mapping from S 2 into R + , verifying S 2 f (x)dx = 1, where dx indicates the Lebesgue measure on S 2 . Let {X i : i = 1, ..., n} be a collection of i.i.d. observations with values in S 2 with common distribution given by f (x)dx. A classical statistical problem, considered for instance by [START_REF] Baldi | Adaptive Density Estimation for directional Data Using Needlets[END_REF][START_REF] Kerkyacharian | Localized Spherical Deconvolution[END_REF][START_REF] Kueh | Locally Adaptive Density Estimation on the Unit Sphere Using Needlets[END_REF], concerns the estimation of f by wavelets/needlets thresholding techniques. To this aim, keeping in mind the notation (3.23), one uses ([10], [START_REF] Hardle | Wavelets, Approximations and Statistical Applications[END_REF]) the following estimator of f :

f (x) = jk β H jk ψ jk (x) , β H jk := β jk I {| β jk | ≥ctn} ,
where t n = log n/n and c is a constant to be determined. Finite-sample approximations on the distributions of β jk can then be instrumental for the exact determination of the thresholding value ct n , see e.g. [START_REF] Donoho | Density estimation by wavelet thresholding[END_REF][START_REF] Hardle | Wavelets, Approximations and Statistical Applications[END_REF].

(Searching for point sources) The joint distribution of the coefficients { β jk } (as defined in (3.23)) is required in statistical procedures devised for the research of so-called point sources, again for instance in an astrophysical context (see for instance [START_REF] Starck | Source detection using a 3D sparse representation: application to the Fermi gamma-ray space telescope[END_REF]). The physical issue can be formalized as follows:

-Under the null hypothesis, we are observing a background of cosmic rays governed by a Poisson measure on the sphere S 2 , with the form of the measure N t (•) defined in (2.6) for some t > 0. In particular, N t is built from a measure N verifying Assumption (2.3), and the intensity of µ t (dx) = E[N t (dx)] is given by the absolutely continuous measure R t • f (x)dx, where R t > 0 and f is a density on the sphere. This situation corresponds, for instance, to the presence of a diffuse background of cosmological emissions.

-Under the alternative hypothesis, the background of cosmic rays is generated by a Poisson point measure of the type:

N * t (A) = N t (A) + P p=1 N (p) t A δ ξp (x)dx ,
where {ξ 1 , ..., ξ P } ⊂ S 2 , each mapping t → N (p) t is an independent Poisson process over [0, ∞) with intensity λ p , and

A δ ξp (x)dx = 1 ⇐⇒ {ξ p ∈ A} .
In this case, one has that N * t is a Poisson measure with atomic intensity

µ * t (A) := E[N * t (A)] = R t A f (x)dx + P p=1 λ p t • A δ ξp (x)dx .
In this context, the informal expression "searching for point sources" can then be translated into "testing for P = 0" or "jointly testing for λ p > 0 at p = 1, ..., P ". The number P and the locations {ξ 1 , ...ξ P } can be in general known or unknown. We refer to [START_REF] Iuppa | A Needlet-Based Approach to the Shower-Mode Data Analysis in the ARGO-YBJ Experiment[END_REF][START_REF] Scodeller | Detection of new point sources in WMAP 7 year data using internal templates and needlets[END_REF] for astrophysical applications of these ideas.

Remark 3.1 In order to directly apply the findings of [START_REF] Peccati | Stein's method and normal approximation of Poisson functionals[END_REF][START_REF] Peccati | Multi-dimensional Gaussian fluctuations on the Poisson space[END_REF], in what follows we shall focus on a randomized version of (3.23), where n is replaced by an independent Poisson number whose parameter diverges to infinity. Also, we will prefer a deterministic normalization over a random one. As formally shown in the discussion to follow, the resulting randomized coefficients can be neatly put into the framework of Section 2.

Needlet coefficients as Wiener-Itô integrals

Let N be a Poisson measure on R + × S 2 satisfying the requirements of Assumption 2.3 (in particular, the intensity of N has the form ρ × ν, where ν(dx) = f (x)dx, for some probability density f on the sphere, and one writes R t = ρ([0, t]), t > 0). For every t > 0, let the Poisson measure N t on S 2 be defined as in (2.6). For every j ≥ 1 and every k = 1, ..., N j , consider the function ψ jk defined in (3.19), and observe that ψ jk is trivially an element of L 3 (S 2 , ν) ∩ L 2 (S 2 , ν) ∩ L 1 (S 2 , ν). We write

σ 2 jk := S 2 ψ 2 jk (x) f (x)dx , b jk := S 2 ψ jk (x) f (x)dx .
Observe that, if f (x) = 1 4π (that is, the uniform density on the sphere), then b jk = 0 for every j > 1. On the other hand, under (2.8),

ζ 1 ψ jk (.) 2 L 2 ≤ σ 2 jk ≤ ζ 2 ψ jk (.) 2 L 2 . (3.24)
Note that (see 3.21) the L 2 -norm of {ψ jk } is uniformly bounded above and below, and therefore the same is true for σ 2 jk (indeed, there exists κ > 0, independent of j and k, such that 0 < κ < ψ jk 2 L 2 (S 2 ) < 1). For every t > 0 and every j, k, we introduce the kernel

h (Rt) jk (x) = ψ jk (x) √ R t σ jk , x ∈ S 2 , (3.25)
and write

β (Rt) jk := Nt h (Rt) jk = S 2 h (Rt) jk (x) Nt (dx) = x∈supp(Nt) h (Rt) jk (x)-R t • S 2 h (Rt)
jk (x)ν(dx) , (3.26) In view of Remark 2.4-(ii), the random variable β (Rt) jk can always be represented in the form

β (Rt) jk = Nt(S 2 ) i=1 ψ jk (X i ) -R t b jk √ R t σ jk ,
where {X i : i ≥ 1} is a sequence of i.i.d. random variables with common distribution ν, and independent of the Poisson random variable Nt (S 2 ). Moreover, the following relations are immediately checked:

E β (Rt) jk = 0, E[( β (Rt) jk ) 2 ] = 1 . (3.27)
Remark 3.2 Using the notation (3.23), we have that

β (Rt) jk = N t (S 2 ) × β (Nt(S 2 )) jk -R t b jk √ R t σ jk .
4

Bounds in dimension one

We are now going to apply the content of Theorem 2.8-( 1) to the random variables β (Rt) jk introduced in the previous section. In the next statement, we write Z ∼ N (0, 1) to indicate a centered Gaussian random variable with unit variance. Recall that

ζ 2 := sup x∈S 2 |f (x)|, p ≥ 1
, and that the constants q p , q ′ p have been defined in (3.21).

Proposition 4.1 For every j, k and every t > 0, one has that

d W β (Rt) jk , Z ≤ (q ′ 3 ) 3 ζ 2 B j √ R t σ 3 jk .
It follows that for any sequence

(j(n), k(n), t(n)), β (Rt(n)) j(n)k(n) converges in distribution to Z, as n → ∞, provided B 2j(n) = o(R t(n) ).
The convergence is σ(N)-stable, in the sense of Proposition 2.10.

Proof. Using (3.25)-(3.26) together with (2.17) and (2.8),

d W β (Rt) jk , Z ≤ S 2 h (Rt) jk (x) 3 µ t (dx) = R t R 3 t σ 3 jk S 2 |ψ jk (x)| 3 f (x) dx ≤ ζ 2 √ R t σ 3 jk ψ jk 3 L 3 (S 2 ) ≤ (q ′ 3 ) 3 ζ 2 B j √ R t σ 3 jk ,
where in the last inequality we use the property (3.21) with p = 3 to have:

ψ jk 3 L 3 (S 2 ) ≤ (q ′ 3 ) 3 B 3j(1-2 3 ) = (q ′ 3 ) 3 B j .
The last part of the statement follows from the fact that the topology induced by the Wasserstein distance (on the class of probability distributions on the real line) is strictly stronger than the topology of convergence in law.

Remark 4.2 For f (x) ≡ {4π} -1 we have

σ 2 jk = 1 4π S 2 ψ 2 jk (x) dx = ψ jk 2 L 2 (S 2 ) ,
and more generally, under (2.8),

d W β (Rt) jk , Z ≤ B j √ R t (q ′ 3 ) 3 ζ 2 ζ 3/2 1 ψ jk 3/2 L 2 (S 2 ) := γ(j, k, t) .
(4.28)

Remark 4.3
The previous result can be given the following heuristic interpretation. The factor B -j can be viewed as the "effective scale" of the wavelet, i.e. it is the radius of the region centred at ξ jk where the value of the wavelet function is not negligible. Because needlets are isotropic, the "effective area" is of order B -2j . For governing measures with density which is bounded and bounded away from zero, the expected number of observations on a spherical cap of radius B -j around ξ jk is hence given by

E card X i : d(X i , ξ jk ) ≤ B -j ≃ R t d(x,ξ jk )≤B -j f (x)dx , ζ 1 B -2j R t ≤ R t d(x,ξ jk )≤B -j f (x)dx ≤ ζ 2 B -2j R t ,
using [4, equation (8)]. Because the Central Limit Theorem can hold only when the effective number of observations grows to infinity, the condition B -2j R t → ∞ is quite expected. In the thresholding literature, coefficients are usually considered up to the frequency J R such that B 2J R ≃ R t / log R t , see for instance [START_REF] Hardle | Wavelets, Approximations and Statistical Applications[END_REF] and [START_REF] Baldi | Adaptive Density Estimation for directional Data Using Needlets[END_REF]; under these circumstances, we have

d 2 β (Rt) J R k , Z = O 1 √ log R t -→ 0 for R t -→ +∞ . Therefore β (Rt)
J R k does converge in law to Z.

Multidimensional bounds

We are now going to apply Part 2 of Theorem (2.8) to the computation of multidimensional Berry-Esseen bounds involving vectors of needlet coefficients of the type (3.26). After having proved some technical estimates in Section 5.1, we will consider two bounds. One is proved in Section 5.2 by means of (2.17), and it is well adapted to the case where the number of needlet coefficients, say d, is fixed. In Section 5.3, we shall focus on (2.16), and deduce a bound which is adapted to the case where the number d is possibly growing to infinity.

A technical result

The following estimate, allowing one to bound the covariance between any two needlet coefficients, will be used throughout this section. We let the notation and assumptions of the previous section prevail.

Lemma 5.1 For any j ≥ 1 and k 1 = k 2 ≤ K j = card{X j } and every τ > 0, there exists a constant C τ > 0, solely depending on τ, and such that

Γ Rt (k 1 , k 2 ) := E β (Rt) jk 1 β (Rt) jk 2 ≤ C τ ζ 2 σ jk 1 σ jk 2 (1 + B j d (ξ jk 1 , ξ jk 2 )) τ .
Proof. We focus on τ > 2; note that the inequality for any fixed value of τ immediately implies the result for all τ ′ < τ. For k 1 = k 2 we have:

|Γ Rt (k 1 , k 2 )| = 1 R t σ jk 1 σ jk 2 S 2 ψ jk 1 (x) ψ jk 2 (x) µ t (dx) = R t R t σ jk 1 σ jk 2 S 2 ψ jk 1 (x) ψ jk 2 (x) f (x)dx ≤ ζ 2 σ jk 1 σ jk 2 S 2 |ψ jk 1 (x)| |ψ jk 2 (x)| dx.
Now we can use a classical argument ( [START_REF] Narcowich | Localized Tight Frames on Spheres[END_REF], [START_REF] Narcowich | Decomposition of Besov and Triebel-Lizorkin Spaces on the Sphere[END_REF][4]) to show that, for any τ > 2, there exists C τ > 0 such that:

|ψ jk 1 |, |ψ jk 2 | L 2 (S 2 ) = S 2 |ψ jk 1 (x) ||ψ jk 2 | (x) dx ≤ κ τ B 2j S 2 1 (1 + B j d (x, ξ jk 1 )) τ 1 (1 + B j d (x, ξ jk 2 )) τ dx .
In order to evaluate this integral, we can for instance follow ( [START_REF] Narcowich | Localized Tight Frames on Spheres[END_REF]), by splitting the sphere S 2 into two regions:

S 1 = x ∈ S 2 : d (x, ξ jk 1 ) > d (ξ jk 1 , ξ jk 2 ) /2 S 2 = x ∈ S 2 : d (x, ξ jk2 ) > d (ξ jk 1 , ξ jk 2 ) /2 .
For what concerns the integral on S 1 , we obtain:

S 1 1 (1 + B j d (x, ξ jk 1 )) τ 1 (1 + B j d (x, ξ jk 2 )) τ dx ≤ 2 τ (1 + B j d (ξ jk 1 , ξ jk 2 )) τ S 1 dx (1 + B j d (x, ξ jk 2 )) τ .
One also has that

S 1 dx (1 + B j d (x, ξ jk 2 )) τ ≤ S 2 dx (1 + B j d (x, ξ jk 2 )) τ = 2π π 0 sin ϑ (1 + B j ϑ) τ dϑ ≤ ≤ 2π B 2j ∞ 0 y (1 + y) τ dy ≤ 2π B 2j 1 0 ydy + ∞ 1 y 1-τ dy ≤ ≤ 2πC B 2j .
Because calculations on the region S 2 are exactly the same and because S 2 ⊂ S 1 ∪ S 2 , we have that, for some constant C τ depending on τ ,

|ψ jk 1 |, |ψ jk 2 | L 2 (S 2 ) ≤ C τ (1 + B j d (ξ jk 1 , ξ jk 2 )) τ ,
yielding the desired conclusion.

Remark 5.2 Assuming that d (ξ jk 1 , ξ jk 2 ) > δ uniformly for all j, we have immediately

E β (Rt) jk 1 β (Rt) jk 2 ≤ κ ′ τ,ζ 2 × B -jτ ,
where the constant κ ′ τ,ζ 2 only depend on τ, ζ 2 .

Remark 5.3

The previous Lemma provides a tight bound, of some independent interest, on the high frequency behaviour of covariances among wavelet coefficients for Poisson random fields. For Gaussian isotropic random fields, analogous results were provided by [START_REF] Baldi | Asymptotics for Spherical Needlets[END_REF], in the case of standard needlets (bounded support), and by [START_REF] Lan | On the Dependence Structure of Wavelet Coefficients for Spherical Random Fields[END_REF]- [START_REF] Mayeli | Asymptotic Uncorrelation for Mexican Needlets[END_REF], in the "Mexican" case where support may be unbounded in multipole space. It should be noted how asymptotic uncorrelation holds in much greater generality for Poisson random fields than for Gaussian field: indeed in the latter case a regular variation condition had to be imposed on the tail behaviour of the angular power spectrum, and in the Mexican case this condition had to be strengthened imposing an upper bound on the decay of the spectrum itself. The reason for such discrepancy is easily understood: for Poisson random fields, non overlapping regions are independent, whence (heuristically) localization in pixel space is sufficient to ensure asymptotic uncorrelation; on the contrary, in the Gaussian isotropic case different regions of the field are correlated at any angular distance, and asymptotic uncorrelation for the coefficients requires a much more delicate cancellation argument.

Fixed dimension

Fix d ≥ 2 and j ≥ 1, consider a fixed number of sampling points {ξ jk 1 , ..., ξ jk d }, and define the associated d-dimensional vector

β (Rt) j• := β (Rt) jk 1 , ..., β (Rt) jk d ,
whose covariance matrix will be denoted by Γ t (note that, by construction, Γ t (i, i) = 1 for every i = 1, ..., d). Our aim is to apply the rough bound (2.17) in order to estimate the distance between the law of β

C -1 op = C 1 2 op = 1 , C -Γ t H.S. ≤ d k 1 =k 2 =1 E β (Rt) jk 1 β (Rt) jk 2 2 ≤ d sup k 1 =k 2 =1,...,d 1 σ jk 1 σ jk 2 C τ ζ 2 (1 + B j d (ξ jk 1 , ξ jk 2 )) τ ≤ d ζ 1 q 2 2 × C τ ζ 2 (1 + B j inf k 1 =k 2 =1,...,d d (ξ jk 1 , ξ jk 2 )) τ = A(t) .
(5.29)

Applying (2.17) yields therefore that

d 2 β (Rt) j , Z ≤ A(t) + d 2 √ 2π 8 d k=1 R t S 2 h (Rt) jk (x) 3 f (x) dx = A(t) + d 2 √ 2π 8
ζ 2 R t R 3 t d k=1 S 2 |ψ jk (x)| 3 σ 3 jk dx ≤ A(t) + d 3 ζ 2 √ R t ζ 3/2 1 q 3 2 √ 2π 8 ψ jk 3 L 3 (S 2 ) ≤ A(t) + (q ′ 3 ) 3 d 3 ζ 2 √ R t ζ 3/2 1 q 3 2 √ 2π 8 B j ,
where we used (3.21) and (3.24) to yield σ 3 jk ≥ ζ 3/2 1 q 3 2 . We write this result as a separate statement.

Proposition 5.4 Under the above notation and assumptions,

d 2 β (Rt) j , Z ≤ d C τ ζ 2 B -jτ ζ 1 q 2 2 (1 + inf k 1 =k 2 =1,...,d d (ξ jk 1 , ξ jk 2 )) τ + (q ′ 3 ) 3 d 3 ζ 2 √ R t ζ 3/2 1 q 3 2 √ 2π 8 B j .
Because τ can be chosen arbitrarily large, it is immediately seen that the leading term in the d 2 distance is decaying with the same rate as in the univariate case, e.g. B j / √ R t . Assuming however that d = d t , i.e. the case where the number of coefficients is itself growing with t, the previous bound may become too large to be applicable. We shall hence try to establish a tighter bound, as detailed in the next section.

Growing dimension

In this section we allow for a growing number of coefficients to be evaluated simultaneously, and investigate the bounds that can be obtained under these circumstances. More precisely, we are now focussing on

β (Rt) j(t)• := β (Rt) j(t)k 1 , ..., β (Rt) j(t)k d t
, where d t → ∞, as t → ∞. Throughout the sequel, we shall assume that the points at which these coefficients are evaluated satisfy the condition:

inf k 1 =k 2 =1,...,dt d ξ j(t)k 1 , ξ j(t)k 2 ≈ 1 √ d t .
(5.30) Condition (5.30) is rather minimal; in fact, the cubature points for a standard needlet/wavelet construction can be taken to form a maximal (d t ) -1/2 -net (see [START_REF] Baldi | Subsampling Needlet Coefficients on the Sphere[END_REF][START_REF] Geller | Continuous Wavelets on Manifolds[END_REF][START_REF] Narcowich | Localized Tight Frames on Spheres[END_REF][START_REF] Pesenson | Sampling of band-limited vectors[END_REF] for more details and discussion). The following result is the main achievement of the paper.

Theorem 5.5 Let the previous assumptions and notation prevail. Then for all τ = 2, 3..., there exist positive constants c and c ′ , (depending on τ, ζ 1 , ζ 2 but not from t, j(t), d(t)) such that we have

d 2 β (Rt) j(t). , Z ≤ cd t 1 + B j(t) inf k 1 =k 2 =1,...,dt d ξ j(t)k 1 , ξ j(t)k 2 τ + √ 2π 8 c ′ d t B j(t) ζ 3/2 1 q 3 2 √ R t . (5.31)
Proof. In view of (2.16) and (5.29), we just have to prove that the quantity √ 2π 8

R t ζ 3/2 1 q 3 2 R 3 t dt k 1 k 2 k 3 S 2 ψ j(t)k 1 (z) ψ j(t)k 2 (z) ψ j(t)k 3 (z) f (z)dz
is smaller than the second summand on the RHS of (5.31). Now note that

dt k 1 k 2 k 3 S 2 ψ j(t)k 1 (z) ψ j(t)k 2 (z) ψ j(t)k 3 (z) dz ≤ λ B(ξ j(t)λ ,B -j(t) ) dt k ψ j(t)k (z) 3 dz ,
where, for any z ∈ B(ξ

j(t)λ , B -j(t) ) dt k ψ j(t)k (z) ≤ dt k C τ B j(t) 1 + B j(t) d(ξ j(t)k , z) τ ≤ C τ B j(t) + dt k:ξ j(t)k / ∈B(ξ j(t)λ ,B -j(t) ) C τ B j(t) 1 + B j(t) d(ξ j(t)k , ξ j(t)λ ) -d(z, ξ j(t)λ ) τ ≤ C τ B j(t) + dt k:ξ j(t)k / ∈B(ξ j(t)λ ,B -j(t) ) C τ B j(t) B j(t) d(ξ j(t)k , ξ j(t)λ ) τ .
Now for ξ j(t)k / ∈ B(ξ j(t)λ , B -j(t) ), x ∈ B(ξ j(t)k , B -j(t) ), we have by triangle inequality

d(ξ j(t)k , ξ j(t)λ ) + d(ξ j(t)k , x) ≥ d(ξ j(t)λ , x),
and because d(ξ j(t)k , ξ j(t)λ ) ≥ d(ξ j(t)k , x), and 2d(ξ j(t)k , ξ j(t)λ ) ≥ d(ξ j(t)λ , x) , we obtain

dt k:ξ j(t)k / ∈B(ξ j(t)λ ,B -j(t) ) C τ B j(t) B j(t) d(ξ j(t)k , ξ j(t)λ ) τ = dt k:ξ j(t)k / ∈B(ξ j(t)λ ,B -j(t) ) 1 meas(B(ξ j(t)k , B -j(t) )) B(ξ j(t)k ,B -j(t) ) κ τ B j(t) B j(t) d(ξ j(t)k , ξ j(t)λ ) τ dx ≤ dt k:ξ j(t)k / ∈B(ξ j(t)λ ,B -j(t) ) 1 meas(B(ξ j(t)k , B -j(t) )) B(ξ j(t)k ,B -j(t) ) κ τ 2 τ B j(t) B j(t) d(ξ j(t)λ , x) τ dx ≤ κ ′ τ B j(t) ,
arguing as in [START_REF] Baldi | Asymptotics for Spherical Needlets[END_REF], Lemma 6. Hence

dt k ψ j(t)k (z) ≤ κ ′′ τ B j(t) , (5.32)
uniformly over z ∈ S 2 , which immediately provides the bound.

λ B(ξ j(t)λ ,B -j ) dt k ψ j(t)k (z) 3 dz ≤ (κ ′′ τ B j ) 3 λ B(ξ j(t)λ ,B -j(t) ) dz = (κ ′′′ B j(t) ) 3 .
Finally, to establish the sharper constraint

S 2 dt k ψ j(t)k (z) 3 dz ≤ κ τ d t B j(t) ,
it is sufficient to note that, exploiting (5.32)

k 1 S 2 ψ j(t)k 1 (z) k 2 ψ j(t)k 2 (z) k 3 ψ j(t)k 3 (z) dz ≤ κ 2 B 2j(t) dt k 1 S 2 ψ j(t)k 1 (z) dz = κ 2 B 2j(t) d j(t) ψ j(t)k L 1 (S 2 ) ≤ d t κ 2 B 2j(t) B -j(t) = d t κ 2 B j(t) ,
where we have used again

ψ j(t)k p L p (S 2 ) = O(B 2j(t)( 1 2 -1 p )p ) = O(B j(t)(p-2)
), for p = 1. Thus (5.31) is established.

For definiteness, we shall also impose tighter conditions on the rate of growth of d t , B j(t) with respect to R t , so that we can obtain a much more explicit bound, as follows:

Corollary 5.6 Let the previous assumptions and notation prevail, and assume moreover that there exists α, β such that, as t → ∞ B 2j(t) ≈ R α t , 0 < α < ), so that even at multipoles in the order of B j(t) = O(10 3 ) we might take around 10 3 coefficients with the multivariate Gaussian approximation still holding. These arrays would not be sufficient for the map reconstruction at this scale, but would indeed provide a basis for joint multiple testing procedures as those described earlier.

Remark 5.8 Assume that d t scales as B 2j(t) ; loosely speaking, this corresponds to the situation when one focusses on the whole set of coefficients corresponding to scale j, so that exact reconstruction for bandlimited functions with l = O(B j ) is feasible. Under this requirement, however, the "covariance" term A(t), i.e. the first element on the right-hand side of (5.31), is no longer asymptotically negligible and the approximation with Gaussian independent variables cannot be expected to hold. The approximation may however be implemented in terms of a Gaussian vector with dependent components. For the second term, convergence to zero when d j(t) ≈ B 2j(t) requires B 3j(t) = o( √ R t ). In terms of astrophysical applications, for R t ≃ 10 12 this implies that one can focus on scales until 180 • /B j ≃ 180 • /10 2 ≃ 2 • ; this is close to the resolution level considered for ground-based Cosmic Rays experiments such as ARGO-YBJ ( [START_REF] Iuppa | A Needlet-Based Approach to the Shower-Mode Data Analysis in the ARGO-YBJ Experiment[END_REF]). Of course, this value is much lower than the factor B j = o( √ R t ) = o(10 6 ) required for the Gaussian approximation to hold in the one-dimensional case (e.g., on a univariate sequence of coefficients, for instance corresponding to a single location on the sphere).

Remark 5.9 As mentioned in the introduction, in this paper we decided to focus on a specific framework (spherical Poisson fields), which we believe of interest from the theoretical and the applied point of view. It is readily verified, however, how our results continue to hold with trivial modifications in a much greater span of circumstances, indeed in some cases with simpler proofs. Assume for instance we observe a sample of i.i.d. random variables {X t } , with probability density function f (.) which is bounded and has support in [a, b] ⊂ R. Consider the kernel estimates

f n (x nk ) := 1 nB -j n t=1 K( X t -x nk B -j ) , (5.34)
where K(.) denotes a compactly supported and bounded kernel satisfying standard regularity conditions, and for each j the evaluation points (x n0 , ..., x nB j ) form a B -j -net; for instance Considering that Nt Rt → a.s. 1, a bound analogous to 5.33 can be established with little efforts for the vector f n (x n. ) := f n (x n1 ), ..., f n (x nB j ) . We leave this and related developments for further research.

  Hess g(x) op , where • op indicates the operator norm. Definition 2.5 The Wasserstein distance d W , between the laws of two random vectors

  of a random Gaussian vector Z ∼ N d (0, I d ), where C = I d stands for the identity d × d matrix. Using Lemma 5.1, one has the following basic estimates:

a

  = x n0 < x n1 ... < x nB j = b , x nk = a + k ba B j , k = 0, 1, ..., B j .As argued earlier, conditionally on N t ([a, b]) = n, (5.34) has the same distribution asf Nt (x nk ) := 1 N t [a, b]B -j b a K( ux nk B -j )dN t (u) ,where N t is a Poisson measure governed by R t × A f (x)dx for all A ⊂ [a, b].

  1 , d t ≈ R β t , 0 < β < α .There exists a constant κ (depending on ζ 1 , ζ 2 , but not on j, d j , B) such that(1 + B j(t) inf k 1 =k 2 =1,...,dt d (ξ jk 1 , ξ jk 2 )) From(5.33), it follows that for R t ≃ 10 12 we can establish asymptotic joint Gaussianity for all sequences of coefficients β

	d 2 β j(t). , Z ≤ κ (Rt)	d t B j(t) √ R t	,		(5.33)
	for all vectors β jk 1 , ..., β (Rt) jk d t (Rt)	, such that (5.30) holds.
	Proof. It suffices to note that
		d t κ ′ τ,ζ 2				1+τ /2 t	)
						= O	d t B j(t) R t √	B (τ +1)2j(t) t R t d τ	1/2
	and
	R t d τ t B (τ +1)2j(t) =	R 1+βτ t R (τ +1)α t	= R t -α+τ (β-α)+1	= o(1) , for τ >	1 -α α -β	.
	Remark 5.7

τ = O(B -τ j(t) d (Rt) j(t)k 1 (t) , ..., β (Rt) j(t)k d (t) ) of dimensions such that d t B j(t) √ R t = o(

1) , e.g. we can take d t ≃ o( √ R t /B j(t) ) ≃ o(10 6 /B j(t)

Research supported by the ERC Grant Pascal, n.277742


Zdeněk Bonaventura

Max Duarte

Anne Bourdon

Marc Massot

Derivation of a merging condition for two interacting streamers in air

The simulation of the interaction of two simultaneously propagating air streamers of the same polarity is presented. A parametric study has been carried out using an accurate numerical method which ensures a time-space error control of the solution. For initial separation of both streamers smaller or comparable to the longest characteristic absorption length of photoionization in air, we have found that the streamers tend to merge at the moment when the ratio between their characteristic width and their mutual distance reaches a value of about 0.35 for positive streamers, and 0.4 for negative ones. Moreover it is demonstrated that these ratios are practically independent of the applied electric field, the initial seed configuration, and the pressure.

Electrical breakdown in air gaps often involves the development of fast ionizing waves that take the form of thin filaments called streamers [START_REF] Raizer | Gas Discharge Physics[END_REF]. These complex and highly nonlinear phenomena precede the formation of sparks and leaders. Similar filamentary structures are observed in sprites: large-scale discharges which appear at low pressure in the altitude range of about 40 to 90 km above large thunderstorms [START_REF] Sentman | [END_REF]3]. Streamers often form branches. Then the individual streamer heads carry charges of the same polarity, and thus they are electrostatically repelling each other. However, careful observations show also the opposite behavior: streamers may even merge or reconnect [4,5,6]. So far, only few numerical studies have been carried out on streamer interactions. First, the impact of charges of numerous positive streamers, propagating simultaneously, on the electric field and on the velocity of one single positive streamer was studied in [7] for a 2D configuration. Later on it was shown in [8] that the electrodynamics of a 2D array of negative streamers substantially differs from the one of a single streamer due to their electrostatic interaction. Recently, the interaction of two streamer discharges in air and other oxygen-nitrogen mixtures has been studied for a 3D cylindrical configuration in [9]. These authors showed that two competing mechanisms have to be considered: attraction due to nonlocal photoionization between streamers and electrostatic repulsion of streamer heads due to their space charges of same polarity. Qualitatively, they concluded that streamer merging is favored if the distance between streamers is smaller or comparable to the longest absorption length of photoionization. However, a more precise or even quantitative description of streamer merging remains an open problem.

In this work we carry out a parametric analysis of the dynamics of interacting streamers to better understand the conditions of streamer merging. Therefore we simulate the interaction of two streamers of the same polarity by means of an accurate numerical method which ensures a time-space error control of the solution. Taking into account that fully 3D simulations are still computationally expensive to carry out an extensive parametric analysis, we have considered a Cartesian 2D geometry to reduce computational costs. This study is then a first step towards a quantitative analysis of streamer merging in real 3D configurations. However, in agreement with [8,10], we consider that the main characteristics of the dynamics of interacting streamers in two dimensions, as described below, will be similar in three dimensions. In practice, we aim at defining geometrical parameters independent of the multivariable complex physical settings and thus a simple characterization that may be extended to more general configurations.

We consider the classical fluid model for air given by drift-diffusion equations self-consistently coupled with Poisson's equation [11,12]:

∂ t n e + ∇ • (n e v e )
+ n n γ + S ph , ∂ t n p + ∇ • (n p v p ) -∇ • (D p ∇n p ) = n e α|v e | -n e n p β ep -n n n p β np + S ph , ∂ t n n + ∇ • (n n v n ) -∇ • (D n ∇n n ) = n e η|v e | -n n n p β np -n n γ,    (1) ε 0 ∇ • E = -q e (n p -n n -n e ), E = -∇V, (2)
where n i is the density of charged species i (e: electrons, p: positive ions, n: negative ions), V and E stand, respectively, for the electric potential and field, and v i = µ i E is the drift velocity. We denote by D i and µ i the diffusion coefficient and the mobility of charged species i, q e is the absolute value of the electron charge, and ε 0 is the permittivity of free space. Moreover α is the impact ionization coefficient, η stands for the electron attachment coefficient, β ep and β np are, respectively, the electron-positive ion and negative-positive ion recombination coefficients, and γ is the detachment coefficient. All these coefficients depend on the local reduced electric field E/N , and thus vary in time and space, where E = |E| is the electric field magnitude, and N is the air neutral density. Further details are given in [15]. For S ph , the photoionization source term for air, we use the 3-Group SP 3 model derived in [13] with Larsen's boundary conditions [14]. In this model, three equivalent absorption lengths are used to model the wavelength dependence of the photoionization source term in air. At ground pressure these absorption lengths are (λ 1 p O2 ) -1 = 0.1408 cm, (λ 2 p O2 ) -1 = 0.0561 cm, and (λ 3 p O2 ) -1 = 0.0105 cm, where p O2 = 158.9 Torr is the partial pressure of molecular oxygen in air.

To handle the stiff modeling equations, we have recently developed in [15] a numerical scheme for multi-scale streamer discharge simulations for general Cartesian multi-dimensional geometries. The strategy is based on a second order time adaptive integration with a splitting technique and dedicated solvers, and space finite volume multiresolution for dynamic grid adaptation. These features involve important efficiency gains in terms of CPU time and memory space while ensuring a time-space error control of the solution. Poisson's equation and Helmholtz equations issued by the photoionization model are considered on the adapted grid, and the resulting sparse linear systems are solved with MUMPS [16,17]. In all simulations carried out in this work (except for some computations in figure 6 where a finer grid is used) a time-space accuracy tolerance of 10 -4 is used for a space resolution of 3.9 µm. This choice guarantees a sufficiently fine time-space representation of the physics, and numerical results disclosing practically the same physical behavior for higher spatial resolutions and tighter accuracy tolerances. The domain size is carefully chosen for the various configurations and with the same space resolution of 3.9 µm such that no interference of the boundaries is evidenced.

In order to preserve similarity constraints at different pressures (corresponding to discharges at different altitudes), all variables are scaled by an appropriate power of the ratio of the air density at ground pressure N 0 , and the air density N at the given pressure [18]. For example, the space dimensions scale as x = x 0 N 0 /N , the charge densities as n = n 0 N 2 /N 2 0 , the electric field as E = E 0 N/N 0 , and the time as t = t 0 N 0 /N , where the subscript "0" refers to the values at ground pressure. Our reference (ground) density for air is N 0 = 2.688 × 10 19 cm -3 at a temperature of 273 K.

We study two separate discharge filaments initiated by placing two identical Gaussian plasma clouds into a homogeneous background electric field generated by two remote electrodes. The pressure, streamer polarity, applied electric field, and the configuration parameters of the initial Gaussian seeds (the width σ, the maximum density n max , and the mutual separation y 0 ) are varied in order to derive a condition for the merging of both streamers.

Figures 1(a)-(c) show snapshots of different phases of the evolution of the net charge density (ρ = n p -n e -n n ) for two simultaneously propagating positive discharge filaments at ground pressure (N/N 0 = 1.0) for an initial seed separation of y 0 = 0.2 N 0 /N cm, a maximum density of n max = 10 13 N 2 /N 2 0 cm -3 , a width σ = 0.02N 0 /N cm, and an applied electric field E a = 1.5E bd , where E bd = 32 N/N 0 kV cm -1 is the breakdown electric field in air. The complete trajectory of the maximum electric field is also shown. Following this trajectory, the transition time t tr , between the streamer repulsion and the beginning of their merging, is defined as the instant when the separation of the two trajectories starts to decrease. For the condition of Figure 1, we have t tr = 6.3 ns. Figure 2(a) shows at this time the distribution of the net charge density. We define in Figure 2(b) the mutual distance of both filaments δ by the separation of the maximum net density ρ max in each head; and the streamer head characteristic width D as the planar length set by the cross section in the direction perpendicular to the applied electric field E a , where the charge density is locally higher than half the maximum net density: ρ ≥ ρ max /2.

We have found that the merging condition for two streamers in air could be represented by the value of the ratio D/δ which appears to be very stable for a given polarity of discharges. To illustrate this, Figure 3 shows the time evolution of the ratio D/δ for two positive or two negative streamers. In this figure, initial separation of the seeds y 0 = 0.15 N 0 /N cm and 0.2 N 0 /N cm are considered, and other parameters are similar to Figure 1. First, Figure 3 shows the time evolution of the ratio D/δ for two negative streamers as these discharges may propagate without photoionization source term. When photoionization is included, the ratio D/δ increases and reaches a value of 0.4 at the moment of transition (4 ns in case of y 0 = 0.15 N 0 /N cm, and 5 ns in case of y 0 = 0.20 N 0 /N cm) between streamer repulsion and the beginning of their merging, and streamer merging is then observed. It is interesting to note that when photoionization is omitted, the ratio D/δ is less than 0.4, and no merging is evidenced. Conversely photoionization must be included for positive streamers as it is necessary for their propagation. In this case Figure 3 shows that the ratio D/δ increases and reaches a value of 0.35 at the transition time 5.6 ns for y 0 = 0.15 N 0 /N cm, (6.3 ns for y 0 = 0.20 N 0 /N cm), and then streamer merging is observed. These results clearly show the crucial role of the photoionization source term on streamer merging. In the following we show that the value of the ratio D/δ at the transition time depends on the polarity of discharges but is practically independent on the multivariable physical settings. For instance, for the same parameters of the seeds as in Figure 1 (n max = 10 13 N 2 /N 2 0 cm -3 , y 0 = 0.2 N 0 /N cm, and σ = 0.02 N 0 /N cm) we have performed a set of simulations with an applied electric field E a varying from 1.5E bd to 2.5E bd at ground pressure. Figure 4 shows that the resulting values of D/δ at the transition times for merging streamers are roughly constant and of about 0.35 for positive polarity discharges, and approximatively 0.4 for negative ones. Hence D/δ is independent of the applied electric field E a for both streamer polarities. Similarly, Figure 5 shows the dependence of D/δ on n max , the maximum density of seeds, for a fixed applied electric field E a = 1.5E bd . Same previous ratios are obtained whereas D/δ is slightly decreasing as n max decreases for positive streamers. However a variation of four orders of magnitude for n max involves only a variation of about 15% of D/δ. In the same way we have tested the variation of D/δ with altitude (not shown here). We have considered altitudes ranging from the ground up to 80 km, corresponding to pressure range from 1 to 1.5 × 10 -5 atm. We have found that the ratio D/δ is independent of the altitude and then of pressure in the range studied.

Moreover we have performed simulations to study the influence of the initial seed separation y 0 on D/δ. Figure 6 shows the dependence of D/δ for y 0 ≤ 0.2 N 0 /N cm. Two sets of seed parameters were considered: (a) for y 0 > 0.1 N 0 /N cm: n max = 10 13 N 2 /N 2 0 cm, σ = 0.02 N 0 /N cm, and E = 1.5E bd ; (b) for y 0 < 0.1 N 0 /N cm: n max = 10 10 N 2 /N 2 0 cm, σ = 0.00736 N 0 /N cm, and E = 2.5E bd . The first set corresponds to the reference condition used for Figure 1. The second set was chosen to be close to the conditions used in [9], and is characterized by a high value of the applied electric field, comparable to typical values of the electric field in streamer heads, and small initial seed width and separation. Because of the small seed width of the second set, for this set, computations were performed for a space resolution of 1.95 µm. Figure 6 shows that for positive streamers with y 0 ≥ 0.04 N 0 /N cm, the value of the ratio D/δ for the second set is slightly less than the reference value of 0.35 for the first set, which can be explained by the small dependence of D/δ on n max , as seen in Figure 5. For negative streamers, Figure 6 shows that for y 0 ≥ 0.04 N 0 /N cm, the ratio D/δ is independent of the value of y 0 . For initial separation of streamers smaller than y 0 < 0.04 N 0 /N cm the transition between repulsion and merging is occurring very close to the initial Gaussian seeds and therefore the region where ρ ≥ ρ max /2 is indistinguishably joined with the charge lobes, typically observed near the initial germs (shown in Figure 2). In this case, the streamer head appears to be significantly wider, and consequently the ratio D/δ increases as shown in Figure 6. It is interesting to note that in the two sets used for Figure 6, very different values of σ have been used. We have carried out additional simulations to verify that the value of σ has a negligible influence on the ratio D/δ.

Finally, for simulations presented in this work, we have found that a stable value for the ratio D/δ is obtained for merging streamers for an initial separation y 0 which is smaller or comparable to the longest absorption length of photoionization for air (0.1408N 0 /N cm). These results clearly show the significant role of photoionization on streamer merging. We have carried out simulations for larger values of y 0 , but it was not possible to derive simple conclusions on streamer interactions from the analysis of the value of the ratio D/δ.

In conclusion, in this work, we have carried out a parametric study on the interaction of two streamers in air using 2D numerical simulations based on an accurate numerical method which ensures a time-space error control of the solution. This study is a first step towards a quantitative analysis of streamer merging in real 3D conditions. We have shown that for a initial separation between streamers that is smaller or comparable to the longest characteristic absorption length of photoionization (0.1408N 0 /N cm), streamers will start to merge, i.e., the distance between their trajectories starts to decrease despite their electrostatic repulsion, at the moment when the ratio D/δ of the streamer characteristic width D and the streamer mutual distance δ attains a value of about 0.35 for positive streamers, and 0.4 for negative ones. These ratios are independent of the applied field, the initial seed configuration, and pressure. Moreover presented results clearly illustrate the significant role of photoionization on the streamer merging. Further studies are required to find the key parameters to characterize streamer interactions when y 0 is considerably larger than the longest photoionization absorption length.

  -∇ • (D e ∇n e ) = n e α|v e | -n e η|v e | -n e n p β ep

Figure 1 :

 1 Figure 1: (color online) Time evolution of the net charge density for two positive streamers at ground pressure (N/N 0 = 1.0) for an applied electric field E a = 1.5E bd and two Gaussian seeds with y 0 = 0.2 N 0 /N cm, n max = 10 13 N 2 /N 2 0 cm -3 , and σ = 0.02N 0 /N cm. (a) t = 5.0 ns: well-developed streamers repulsing each other, (b) t = t tr = 6.3 ns: transition between repulsion and merging, (c) t = 8.0 ns: propagation of a single discharge. Black solid line: trajectory of the maximum electric field.

Figure 2 :

 2 Figure 2: (color online) (a) Net charge density and (b) Region where the net charge density is higher than half the maximum net charge density at t = t tr for the same condition as in Figure 1(b). Definition of the streamer head characteristic width D and the mutual distance of both filaments δ.

Figure 3 :Figure 4 :Figure 5 :

 345 Figure 3: (color online) Time evolution of the ratio D/δ for (a) two negative streamers with photoionization source term, (b) two negative streamers without photoionization source term, (c) two positive streamers with photoionization source term. Results for initial separation of the seeds y 0 = 0.15 N 0 /N cm and y 0 = 0.2 N 0 /N cm are shown. Other parameters are similar to Figure 1. Critical merging ratios are indicated with stripes.

Figure 6 :

 6 Figure 6: (color online) Ratios D/δ for positive and negative streamers for different initial seed distances y 0 , with (a) n max = 10 13 N 2 /N 2 0 cm, E = 1.5E bd , σ = 0.02 N 0 /N cm; (b) n max = 10 10 N 2 /N 2 0 cm, E = 2.5E bd , σ = 0.00736 N 0 /N cm. Absorption lengths for photoionization at ground pressure are shown.

Acknowledgements

This research was supported by a fundamental project grant from ANR (French National Research Agency): Séchelles (ANR-09-BLAN-0075-01, PI. S. Descombes), and by a DIGITEO RTRA project: MUSE (PI. M. M.). Z. B. acknowledges support by project CZ.1.05/2.1.00/03.0086 funded by European Regional Development Fund and support of Ecole Centrale Paris.