Présentée Devant

Liang Yu

Ivan Markovsky

Bert Roozen

Manuel Melon

Barbara Nicolas

Professeur) LASPI François Guillet

Roanne Examinateur

Jérôme Antoni

R Gourdon

M Jean

Marc Lancelin

M Gérard Scorletti
email: gerard.scorletti@ec-lyon.fr

H Charles

Gudrun Mme

Bornette

Safia Ait

Chalal Insa

M Lagarde

M Didier

M Johannes Kellendonk

M Laboune

Jean-Yves Buffiere
email: yves.buffiere@insa-lyon.fr

Mega

M Philippe Boisse

Scso Scso

M Obadia

Lionel Sec

Viviane Polsinelli

J Y Toussaint

Michael Vannier

Thibault Lafont

Marion Berton

Sophie Baudin

Roch Scherrer

Tit For Tat Thibaut Le Magueresse

Corentin Chesnais

Valentin Meyer

Youssef Gerges

Fulbert Mbailassem

Meriem Dib

Edouard Cardenas Cabada

Geng Lin

Su Ting

Yao Min

Sandra Raissa

Loic Forget

Souha Grau

Christophe Kassab

Benjamin Marchetto

William Trevisan

Gousseau

Meriem Elazami

Iulia Mirea

Younes Farouj

Claire Mouton

Mohammad Azizian Kalkhoran

Han Lei

Jia Tan

Yin Xun Neng

Qian

Yin Wang

Yin

Bi Xiong

Jing

Xue Jiang

Liu Jiao

Ji Xun

Yan Xi Min

Zhan Bo

Wu Zhao

Acoustical source reconstruction from non-synchronous sequential measurements

Keywords: inverse acoustic problem, sequential measurements without reference, sparse eigenvalue spectrum models, cyclic projection xii spectral matrix completion Methodologie: sparse eigenvalue spectrum models Research Background: sequential measurements Formulation

be considered as a better alternative to the classical referenced methods when the number of references is insufficient or references are of poor quality ; (3) CP is suggested to be used in cases with low SNRs and relatively small numbers of sources, while FISTA is suggested to be used in cases with high SNRs and relatively large numbers of sources.

I write this acknowledgement part sincerely and humbly since I deeply realized that this thesis will not be finished without the assistances and kindness of my mentors, my family, my friends, my colleagues ... First, I would like to warmly thank my supervisors Jerome Antoni and Quentin Leclere. Mind is not a vessel that needs filling but wood that needs igniting, they are exactly the persons that light my mind up. I would like to thank the professors of the jury for my PHD thesis : Ivan Markovsky, Bert Roozen, Manuel Melon, Barbara Nicolas, François Guillet. It is their critical comments that really improve the quality of this PHD thesis. I would like to thank the continuous discussions with Dr. Antonio Pereira and Professor Ivan Markovsky. Dr. Antonio Pereira is not only a good colleague to work together, but also a good friend to share the life. I would like to thank my parents who give me infinite courage and eternal supports to meet every challenge in my life. I would like to thank my master thesis supervisors Professor Fang Yong (ShangHai University, China) and Professor Kamata sei-ichiro (Waseda University, Japan) who encourage me to continue to work in signal processing (I was fascinated and ambitious to understand all the signal processing techniques : Fourier Transform, Hilbert-Huang transform, wavelet, dictionary learning, sparse representation ...). I would like to thank my French teacher Corinne Lotto and her mother for their patience and warmth during my stay in France. I have met really a lot of cool friends in Lyon and it is them that make my life so wonderful : Rainer Stelzer (technique geek and musician), Jeremy Moriot (crazy party guy), Du Liang Fen (Chinese girl No.1 ; No.3 for Anders), Ha Dong Hwang (da ge and he has seen the world), Li Sheng Fu (best comedian), Anders Lindberg (depression), Razvan Stoica (encyclopaedia), Zheng Li Lei (machine learning geek), Lei Xiao Qin (3 years jogging), Li Xin Hui (vibration geek) ... I would like to thank all the professors, secretaries, students, staffs, engineers of LVA, v

Résumé

Une limitation fondamentale du problème inverse acoustique est déterminée par la taille et la densité de l'antenne de microphones. Une solution pour atteindre un grande antenne et / ou à forte densité de microphones est de scanner l'objet d'intérêt en déplaçant séquentiellement une antenne de petites dimensions (mesures dites séquentielles). La différence entre des mesures séquentielle et une mesure simultanée en tous points est que dans ce dernier cas, l'intégralité de la matrice interpectrale peut être estimée, contrairement au cas des mesures séquentielles qui ne permettent l'estimation que d'une partie réduite de cette matrice ; les éléments croisés (interspectres) entre deux points de mesures n'appartenant pas à la même séquence ne sont pas estimés. Néanmoins, ces données restent nécessaires pour la reconstruction acoustique. Dans l'approche classique, une ou plusieurs références sont utilisées pour retrouver les données manquantes. L'objet de cette thèse est de récupérer les éléments manquants de la matrice interspectrale sans capteurs de référence, dans le cas où le champ acoustique est suffisamment cohérent pour mettre en ouvre les mesures séquentielles. Deux modèles de spectre de valeurs propres parcimonieux sont proposés pour résoudre ce problème, le premier impose une valeur faible de rang, tandis que le second repose sur la minimisation de la norme nucléaire (recherche d'une solution faiblement parcimonieuse).

Les modèles proposés sont construits sur deux fondements physiques : le faible rang de matrice spectrale et la continuité du champ acoustique. Le premier modèle, appelé structured low rank model, (modèle de rang réduit sous contrainte) se résume à la recherche d'une matrice interspectrale de rang réduit, conforme aux données mesurées, observant une certaine continuité spatiale du champ de pression, et conservant les propriétés fondamentales de symétrie hermitienne et de positivité. Un algorithme de Projection Cyclique (CP-cyclic projection) est proposé dans ce travail pour trouver une solution optimale à l'intersection entre trois ensembles prédéfinis. Le second modèle, appelé weakly sparse eigenvalue spectrum (spectre de valeurs propres faiblement parcimonieux) se résume à ix minimiser la norme nucléaire de la matrice spectrale soumise aux données mesurées, observant une certaine continuité spatiale du champ de pression, et conservant les propriétés fondamentales de symétrie hermitienne et de positivité. Un algorithme iteratif rapide de seuillage-troncature (Fast iterative shrinkage thresholding algorithm -FISTA) est proposé pour résoudre ce problème de minimisation de norme nucléaire.

Les méthodes proposées sont analysées et comparées par des simulations numériques de diverses configurations et sont également validées expérimentalement. Trois principaux résultats sont obtenus dans cette thèse :

(1) CP et FISTA peuvent être à la fois utilisés dans le cas de mesures séquentielles avec ou sans références ; (2) CP et FISTA peuvent être considérés potentiellement comme une meilleure alternative à l'approche classique (avec références) ; (3) CP est proposé pour être utilisé dans les cas de faible RSB (rapport signal à bruit), et un petit nombre de sources, tandis que FISTA est proposé pour être utilisé dans le cas d'un RSB élevé et un grand nombre de sources.

Mots clés : reconstruction inverse acoustique, mesures séquentielles sans référence, modèles de spectre de valeurs propres parcimonieux, projection cyclique, x

Abstract

A fundamental limitation of the inverse acoustic problem is determined by the size of the array and the microphone density. A solution to achieve large array and/or high microphone density is to scan the object of interest by moving sequentially a small prototype array, which is referred to as sequential measurements. In comparison to a large array and/or high microphone density array that can acquire simultaneously all the information of the spectral matrix, in particular all cross-spectra, sequential measurements can only acquire a block diagonal spectral matrix, while the cross-spectra between the sequential measurements remain unknown due to the missing phase relationships between consecutive positions. Nevertheless, these unknown cross-spectra are necessary for acoustic reconstruction. The object of this thesis is to recover the missing elements of the spectral matrix in the case that the acoustical field is highly coherent so as to implement the sequential measurements. Sparse eigenvalue spectrum are assumed to solve this problem, which lead to a structured low rank model and a weakly sparse eigenvalue spectrum model. The proposed models are constructed upon two physical pillars : the low rank of the spectral matrix and the continuity of the acoustical field. Structured low rank model is shown to boil down to finding a full spectral matrix subject to reduced rank, measurements fitting, spatial continuity of the sound field, and constraint of hermitian symmetry. A Cyclic Projection (CP) algorithm is proposed in this work to find an optimal solution at the intersection between three predefined sets. Weakly sparse eigenvalue spectrum model is shown to boil down to minimizing the nuclear norm of the spectral matrix subject to measurements fitting, hermitian symmetry and spatial continuity of the sound field. A Fast Iterative Shrinkage Thresholding Algorithm (FISTA) is proposed to solve this nuclear norm minimization problem. The proposed methods are analysed and compared through numerical simulations under diverse setups and are also validated experimentally. Three main results are obtained in this thesis : [START_REF] Fahy | Foundations of engineering acoustics[END_REF] CP and FISTA can be both used to realize sequential measurements without and with few references ; Inverse acoustic reconstruction consists in measuring the acoustic quantity in the field by a sensor array (forward problem) and reconstructing the image of noise sources by back-propagation algorithms (inverse problem), which has a wide range of applications in machine fault diagnosis and noise control of mechanical equipments etc. A fundamental limitation of inverse acoustic reconstruction is determined by the size of the array and the microphone density. Specifically, the minimum spatial resolution hinges on the size of the array and the maximum spatial resolution depends on the maximum distance between the microphones. Thus the working frequency range is

c min{L x , L y } ≤ f ≤ c 2max{d x , d y } , (1.1)
where c is the sound speed and L x , L y , d x , d y are shown in Fig. (1.1) (a). For the sake of improving the spatial resolution of inverse acoustic reconstruction, sequential measurements are investigated in this thesis, which is a solution to achieve large array and/or high microphone density by moving sequentially a small prototype array to scan the object of interest. For instance, Fig. 1.1 (a) shows a prototype array that is optimized to a given configuration as compared with the source dimensions, Fig. 1. 1 (b) shows how to use sequential measurements to make the array denser, and Fig. 1.1 (c) how to make the array larger. It is defined that the spectral matrix is a matrix whose element is the covariance between the measured pressures (a vector concatenate either with finite or infinite number of snapshots) at corresponding two measurement positions in an array for a given frequency. In comparison to a large array and/or high microphone density array that can acquire simultaneously all the information of the spectral matrix, in particular all cross-spectra, sequential measurements can only acquire a block diagonal spectral matrix, while the cross-spectra between the sequential measurements remain unknown due to the missing phase relationships between consecutive positions. Nevertheless, these unknown cross-spectra are necessary for acoustic reconstruction. The object of this thesis is to propose methods to recover the missing elements of the spectral matrix (as would be measured by simultaneous measurement). It assumes only the knowledge of the spectral sub-matrices Ŝ(i) pp , i = 1, ..., P with size M × M , where P is number of measurements and M is number of microphones in the array. Let Ŝm pp be the spectral matrix constructed by rearranging all spectral sub-matrices Ŝ(i) pp , i = 1, ..., P methodically in block diagonal positions and by padding the remaining positions with zeros (as shown on the left side of Fig. 1.2 ; the spectral matrix returned by sequential measurements Ŝm pp (with size M P × M P ) is related to the unknown spectral matrix S (with size M P × M P , and which is shown on the right side of Fig. 1.2 as follows

A (S) + E = Ŝm pp , ( 1.2)
where the sampling operator A ∶ C M P ×M P → C M P ×M P extracts the elements from diagonal block matrix parts, and the error E is composed of estimation error and measurement error. Then the problem can be formulated generally as searching for a full spectral matrix S given the partial spectral matrices Ŝm pp , which is a matrix completion problem. The assumption of a sparse eigenvalue spectrum(sparsity implies that only a few components in eigenvalue spectrum are non-zero and the rest are zero) is proposed to solve this problem, which leads to a structured low rank model and a weakly sparse eigenvalue spectrum model. The physical basis underlaying this approach is founded upon two properties : the low rank of the spectral matrix and the continuity of acoustic field.

1. Low rank spectral matrix. The high correlation of the acoustical field implies a low rank of the spectral matrix. The rank of the spectral matrix measures the stochastic model complexity of the acoustic field, which is consistent with the number of "virtual sources"(equivalent uncorrelated sources that compose the acoustical field) required to produce the pressure measurement. For instance, a rank one spectral matrix corresponds to a fully coherent field (the measurement in one position can be easily predicted from other positions in the field), whereas a spectral matrix with infinite rank corresponds to a diffuse field [START_REF] Fahy | Foundations of engineering acoustics[END_REF] in general.

2. Spatial continuity of the acoustic field. The position coordinates of sequential measurements reflect the acoustic field correlation due to the continuity of the acoustic field. For example, when the distance between two consecutive measurements are very close to each other, the sensed fields will be very similar which indicates the high correlation between the sequential measurements. Unfortunately, the information on position coordinates is not explicated in the spectral matrix, which will have to be introduced by other means.

These two specific properties imply the high correlation of the acoustical field which is considered as the fundamental physical assumption to formulate the problem.

Organization of the thesis

The contributions and contents of each chapter of this thesis are explained separately as follows, and a big picture about this thesis is shown in Fig. 1.3. In chapter 2, entitled "Literature survey", a review of the various techniques that have been proposed to solve the inverse acoustic problem is made ; three categories of methods are surveyed in this thesis, namely : (1) Single measurement, (2) Sequential measurements with references, (3) Sequential measurements without reference. In particular, methods based on assumptions of band-limitation or sparsity are explained, which both imply a high correlation of the acoustical filed. The perspectives for each category and their connections with the methods proposed in this thesis are also discussed.

In chapter 3, entitled "Structured low rank model and Cyclic Projection", the problem of recovering the missing elements in the spectral matrix is formulated as a structured low rank model, and is shown to boil down to finding a spectral matrix subject to a given constraint of hermitian symmetry, measurements fitting, reduced rank, and spatial continuity of the acoustic field. A Cyclic Projection (CP) algorithm is proposed correspondingly to find an optimal solution at the intersection between three predefined sets in order to solve this structured low rank problem.

In chapter 4, entitled "Weakly sparse eigenvalue spectrum model and fast iterative shrinkage thresholding algorithm", the eigenvalue spectrum of the spectral matrix is assumed to be weakly sparse which means the spectral matrix is allowed to be of full rank but with only a few dominant eigenvalues. Thus, the nuclear norm is used as a measure of a weakly sparse eigenvalues spectrum (i.e. the sum of eigenvalues). Not only is this assumption more general than the previous one, but it also avoid the difficulty of a priori determining the actual rank of the spectral matrix. Fast iterative shrinkage thresholding algorithm (FISTA) is proposed correspondingly to solve this weakly sparse eigenvalue spectrum problem.

In chapter 5, entitled "Simulation and experimental validation", a common simulation platform is constructed with various simulation setups, which is used to separately test CP and FISTA. The different simulation scenarios are constructed as follows : sequential measurements without reference, sequential measurements with various SNRs, sequential measurements with different shift distances, and sequential measurements with references.

In chapter 6, entitled "Exploring the difference between CP and FISTA", some simulation scenarios are constructed to illustrate the different mechanisms between CP and FISTA with various SNRs and numbers of sources in order to simulate different eigenvalue spectra. The performance analysis is conducted both from the theoretical and practical viewpoints. The differences between the two proposed algorithms are investigated and a guideline to choose the algorithm is given. Acoustical reconstruction aims at characterizing acoustical sources that radiate at certain distance from partial measurements of the pressure field, which is considered as an inverse problem [START_REF] Tarantola | Inverse problem theory and methods for model parameter estimation[END_REF] [START_REF] Santamarina | Discrete signals and inverse problems : an introduction for engineers and scientists[END_REF]. Classical methods consist in interpolating the measurements onto some spatial basis, implicitly or explicitly, and then the source field is reconstructed by back-propagation that relies on the principle of Least Square (LS), Maximum Likelihood (ML), Maximum a Posteriori (MAP). The choice of a spatial basis varies with the topology of the source surface, the a priori spatial distribution of the sources, the geometry of the array, the type of propagation (near-field or far-field, free field or scattering), and the frequency range of interest etc [START_REF] Antoni | A bayesian approach to sound source reconstruction : Optimal basis, regularization, and focusing[END_REF]. Beamforming [START_REF] Johnson | Array signal processing : concepts and techniques[END_REF] assumes that the discrete sound field has a sparse distribution of well-separated monopoles (note that beamforming is generally a localization method rather than a reconstruction method). The Equivalent source methods (ESM) [START_REF] Leclere | Acoustic imaging using under-determined inverse approaches : Frequency limitations and optimal regularization[END_REF] [7] [START_REF] Pereira | A theoretical and experimental comparison of the equivalent source method and a bayesian approach to noise source identification[END_REF] assumes the sources are modeled by a dense distribution of monopoles. The Near-field acoustical holography (NAH) [START_REF] Maynard | Nearfield acoustic holography : I. theory of generalized holography and the development of nah[END_REF] assumes that the sources can be expanded as the Fourier basis. Steiner and Hald [START_REF] Steiner | Near-field acoustical holography without the errors and limitations caused by the use of spatial dft[END_REF] [START_REF] Hald | Basic theory and properties of statistically optimized near-field acoustical holography[END_REF] proposed a statistically optimal near-field acoustical holography (SONAH), which the sound field is expressed as a linear superposition of plane waves and evanescent wave (without using the Fast Fourier Transform (FFT) as in NAH). The Inverse boundary element method (IBEM) [START_REF] Schuhmacher | Sound source reconstruction using inverse boundary element calculations[END_REF] [13] [14] [START_REF] Borgiotti | The power radiated by a vibrating body in an acoustic fluid and its determination from boundary measurements[END_REF] is based on numerical transfer function obtained by the boundary element method (BEM). The Helmholtz equation-least-squares method (HELS) [START_REF] Wang | Helmholtz equation-least-squares method for reconstructing the acoustic pressure field[END_REF] [START_REF] Wu | On reconstruction of acoustic pressure fields using the helmholtz equation least squares method[END_REF] represents the pressure field as an expansion of spherical harmonics. The Wave superposition algorithms (WSA) [START_REF] Li | A sound field visualization system based on the wave superposition algorithm[END_REF] [19] [20] [21] [22] use a simple source (monopole, dipole, etc.) as the equivalent sources to get the source strength density (weighting factor). The spatial bases used by different methods are listed in table 2.1. These classic methods measure the complex sound pressure field that is located near the radiation source, and are able to reconstruct the entire three-dimensional acoustic space (such as sound pressure, particle velocity and sound intensity). When measurements are close to the sound source in the near-field region, they can capture evanescent waves whose amplitudes decay exponentially with the propagation distance and which contains a wealth of sound source details, thus breaking Rayleigh's limit [START_REF] Williams | Holographic imaging without the wavelength resolution limit[END_REF] as encountered in earlier holography techniques [START_REF] Metherell | Introduction to acoustical holography[END_REF] [START_REF] Shewell | Inverse diffraction and a new reciprocity theorem[END_REF]. The upper plot in Fig. 2.2 shows the spatial sampling in NAH : if the array aperture covers the whole acoustical field, then the spatial signal can be perfectly reconstructed according to the Nyquist-Shannon sampling theorem provided that k s ≥ 2k m [START_REF] Johnson | Array signal processing : concepts and techniques[END_REF], where k s is the spatial sampling frequency and k m is the acoustic field bandwidth as shown in the left figure of Fig. 2.3.

Literature Survey

Inverse acoustic problem

Approach Spatial basis

NAH

plane waves/cylindrical harmonics/spherical harmonic (regular spatial sampling) IBEM generalized spherical harmonics (or "acoustical modes" in general) ESM dense distribution of monopoles WSA monopoles and/or dipoles (inside the object) HELS spherical harmonic SONAH plane waves/cylindrical harmonics/spherical harmonic Table 2.1 -Classic methods and their corresponding basis to be used Nevertheless, the measurements in the above methods are taken with a microphone array whose dimensions and density are circumscribed due to obvious cost and technical reasons ; this brings about strong limitations on the reconstruction capability in terms of spatial resolution and magnitude quantification. Take as an example in Near-field acoustical holography (NAH) [START_REF] Williams | Fourier acoustics : sound radiation and nearfield acoustical holography[END_REF] where it is assumed that the measurement aperture is infinite in theory for the purpose of guaranteeing the reconstruction quality ; however, a limited measurement aperture is only available to cover a partial radiating region as shown in the middle figure of Fig. 2.2, which results in truncation effects, wrap-around errors and other errors due to the discontinuity of the aperture edge and finite aperture effects [START_REF] Williams | Fourier acoustics : sound radiation and nearfield acoustical holography[END_REF]. In order to reduce the errors resulting from the finite aperture effects and to ensure the accuracy of the holographic reconstruction, NAH places a strong restriction on the measurements aperture size, such that the holographic aperture should be at least four times as large as the size of the sound source, so that the sound pressure outside of the measurement aperture is allowed to reach a negligible level. This is possible for small size (with respect to the array) sources, however it is unrealistic for sound sources with large dimensions (with respect to the array) such as vehicle, submarine and machinery etc. The second issue is spatial resolution, which is determined by the density of the sensors in the array. Thus, spatial resolution is one of the most important index criteria of acoustic imaging, which determines the amount of information an image can provide.

In order to breakthrough these crucial limitations, different techniques have been proposed to ameliorate the reconstruction quality, which can be categorised in three groups based on the setup of measurements, namely : Single measurement, Sequential measurements with references, Sequential measurements without reference. Our goal in this chapter is to survey the state-of-the-art of acoustic imaging techniques devoted to pushing the aforementioned limits, and their relation with the methods proposed in this thesis as illustrated in Fig. 2.1.

Single measurement

Single measurement is defined as the configuration where all the measurements are captured at once and simultaneously. Two categories of techniques are surveyed : [START_REF] Fahy | Foundations of engineering acoustics[END_REF] Patch near-field acoustical holography (Patch NAH) and (2) Sparsity induced acoustical source reconstruction.

Patch near-field acoustical holography

As a result of the limited aperture, conventional NAH will generally produce reconstruction error which may be detrimental in practice. The Patch NAH employs numerical method to extrapolate the sound pressure data outside the measurements aperture, so as to simulate a larger measurement aperture, and this indirectly increases the measurement aperture size and reduces the finite aperture effect. The measured sound pressure is forced to be continuous at the edge of the aperture which may considerably reduce reconstruction errors. Patch NAH breaks the limitation of the aperture size to some degree and reduces somewhat the amount of measurements. Saijyou and Yoshikawa [START_REF] Saijyou | Reduction methods of the reconstruction error for large-scale implementation of near-field acoustical holography[END_REF] investigated two kinds of extrapolation method, namely k-space data extrapolation and real-space data extrapolation. The k-space data extrapolation is implemented by zeros-padding hologram surface and iteratively filtering in the wave-number domain ; real-space data extrapolation is implemented by zeros-padding the hologram surface and iteratively forward and back propagation the acoustical field. Williams [START_REF] Williams | Continuation of acoustic near-fields[END_REF] reformulated this extrapolation process as an inverse problem and a modified Tikhonov regularization step was developed. It assumes that the coherence of the pressure field enables the extension of the measured field into the region that lies outside the measured aperture, as a "continuation" of the acoustical field. At the same time, a Fast Fourier Transform (FFT) based and singular value decomposition based patch NAH [START_REF] Williams | Fast fourier transform and singular value decomposition formulations for patch nearfield acoustical holography[END_REF] were proposed to deal with planar surfaces and general shapes. Subsequently, Lee and Bolton [START_REF] Lee | Patch near-field acoustical holography in cylindrical geometry[END_REF] generalized the FFT based patch NAH from planar surfaces to cylindrical geometry. It was found that the zero-padding processes in cylinder patch NAH and planar patch NAH are different since a cylindrical sound source is circumferentially closed and thus the circumferential pressure field should also be periodic. The FFT based patch NAH has an advantage of simplicity and rapidity in extrapolation algorithms, nevertheless it is only applicable to the case of uniform and regular surfaces of orthogonal coordinate system, thus Saijyou and Uchida [START_REF] Saijyou | Data extrapolation method for boundary element method-based near-field acoustical holography[END_REF] proposed a boundary element method (BEM) based patch NAH which is applicable to the objects with any structured shape. Correspondingly, Sarkissian [START_REF] Sarkissian | Method of superposition applied to patch near-field acoustic holography[END_REF] proposed the equivalent source methods (ESM) based patch NAH and Pascal et al. [START_REF] Pascal | Patch near-field acoustic holography : Regularized extension and statistically optimized methods[END_REF] introduced the statistically optimized methods based patch NAH. Zhang et al. [START_REF] Zhang | Patch near field acoustic holography based on particle velocity measurements[END_REF] investigated patch NAH based on particle velocity measurements and claimed that it can achieve better performance since the particle velocity decays faster toward the edges of the measurements aperture than pressure. Lee and Bolton [START_REF] Lee | A one-step patch near-field acoustical holography procedure[END_REF] proposed an alternating orthogonal projection algorithm implemented as a "smooth-and-replace" two step iteration. The theoretic analysis of the algorithm was given in [START_REF] Lee | Reconstruction of source distributions from sound pressures measured over discontinuous regions : Multipatch holography and interpolation[END_REF] and a multi-patch holography procedure was proposed correspondingly to generalize the previous work. Scholte et al. [START_REF] Scholte | Truncated aperture extrapolation for fourier-based near-field acoustic holography by means of border-padding[END_REF] extrapolated the pressure fields by a linear predictive filter using the measured pressures ; then a linear predictive border-padding method was used as an extrapolation technique that greatly reduced the leakage and spatial truncation errors in planar NAH. Reference [START_REF] Bissinger | Violin f-hole contribution to far-field radiation via patch near-field acoustical holography[END_REF] investigated the violin f-hole contribution of the far-field radiation by FFT based patch NAH, and it is proved that patch NAH has a good local characterization ability which can separate the f-hole contribution from the surface motion contributions to the total radiation of the corpus below 2.6kHz. Discussion : Patch NAH consists in extrapolating the acoustic field outside the aper-ture from the hologram data that is measured within the aperture relying on the bandlimited assumption of the spatial spectrum. In sparsity induced acoustical source reconstruction, the acoustical source is modeled as a combination of a few elements taken in an over-completed dictionary ; the assumption of sparsity permits to have a more compact representation than usually allowed by a predefined basis. Patch NAH can be recognized as the application of band-limited extrapolation methods in acoustic. This is a classic signal processing problem which has been addressed since a long time, connected with reconstructing the whole signal from only partial signal measurements [START_REF] Chen | Some aspects of band-limited extrapolations[END_REF] [START_REF] Sanz | Some aspects of band-limited signal extrapolation : Models, discrete approximations, and noise[END_REF]. Popular solutions are given by the Papoulis-Gerchberg algorithm [START_REF] Papoulis | A new algorithm in spectral analysis and band-limited extrapolation[END_REF] and its enhanced versions [START_REF] Salari | Super-resolution using an enhanced papoulis-gerchberg algorithm[END_REF]. A basic assumption of patch NAH is that the acoustic field is band limited in the wave-number space (K-space) as shown by the red curve in the middle of the Fig. 2.3. This assumption implies the high coherence of the acoustical field. In the practical case where the measured signal has a non-finite spatial spectrum as shown by the black curve in the middle of Fig.

2.3, Patch NAH can only recover a low pass version of actual signal.

Sparsity induced acoustical source reconstruction

In Ref [START_REF] Bai | Application of convex optimization to acoustical array signal processing[END_REF], the sparse assumption for the spectrum is investigated instead of the bandlimited assumption, which provides another perspective for extrapolating the acoustical field. The acoustical field was firstly approximated by an expansion on a few plane waves, then the field interpolation problem was formulated as a constrained l 1 norm optimization which was solved by CVX (it is a Matlab-based modeling system for convex optimization). Furthermore, a robust Bayesian super-resolution approach was proposed by Chu et al. [START_REF] Chu | Robust bayesian super-resolution approach via sparsity enforcing a priori for near-field aeroacoustic source imaging[END_REF] : an improved forward model of aeroacoustic power propagation was introduced, which accounts for both background noise and forward model uncertainty ; then a sparse prior for source powers was enforced by means of a double exponential distribution. The Bayesian framework allows both the hyper-parameters and the source powers to be estimated by the joint Maximum A Priori (MAP). Near-field acoustical holography using sparse regularization and compressive sampling principles was proposed by Chardon et al. [START_REF] Chardon | Nearfield acoustic holography using sparse regularization and compressive sampling principles[END_REF]. An over-complete dictionary was built for convex homogeneous plates with arbitrary boundary conditions, then the source normal velocities was decomposed as the sum of a few atoms in the dictionary ; by using the compressed sensing framework, the sparse signals could be reconstructed with significantly less measurements than conventional methods. Suzuki [46] [47] proposed the Generalized Inverse Beamforming (GIB) that aims at identifying coherent or incoherent, distributed or compact, monopole or multipole sources. The source identification problem was formulated as a l 1 norm minimization problem where the source distribution is modeled as a linear combination of a series of monopoles and dipoles. The GIB method was also used in the test of air jets [START_REF] Dougherty | Improved generalized inverse beamforming for jet noise[END_REF], where rectangular jets can have asymmetric radiation patterns. In the context of acous-tic measurements of enclosed spaces, Pereira [START_REF] Pereira | Acoustic imaging in enclosed spaces[END_REF] [START_REF] Oudompheng | A theoretical and experimental comparison of the iterative equivalent source method and the generalized inverse beamforming[END_REF] proposed an iterative Equivalent Source Method (iESM), based on iterative weighted least squares, which also promotes the sparsity of the source in some sense. Malioutov et al. [START_REF] Malioutov | A sparse signal reconstruction perspective for source localization with sensor arrays[END_REF] proposed a source localization method based on sparsity-enforcing regularization ; the approach involves reformulating the source location problem in a constrained optimization by using an overcomplete basis and sparsifying regularization, which results in focusing the signal energy an an increase of the resolution. Co-sparse signal model (analysis model) can be seen as an analytical version of the classic sparse-land model (synthesis model) [START_REF] Elad | Sparse and redundant representations : from theory to applications in signal and image processing[END_REF], and it has also been used in sound source localization [START_REF] Kitić | Hearing behind walls : localizing sources in the room next door with cosparsity[END_REF]. Kitic [START_REF] Kitić | A review of cosparse signal recovery methods applied to sound source localization[END_REF] has compared several state-of-the-art methods for co-sparse signal recovery in the context of sound source localization.

Discussion : In sparsity induced acoustical source reconstruction, it is assumed that the measured spatial acoustic signal can be expanded with a few atoms in an over-completed dictionary as shown in the right part of Fig. 2.3, A sparse representation also implies a high correlation of the acoustical field, yet it is more general than the band-imitated assumption which is a particular case. Thus it provides a powerful approach to address the reconstruction of a acoustical sources with limited (missing) measurements.

Sequential measurements with references

A solution to achieve large array and/or high microphone density is to scan the object of interest by moving sequentially a small prototype array, which is referred to as sequential measurements. A small prototype array is moved along different locations in front of or around the radiating object to augment the microphone density or the measurement aperture size. For instance, Fig. 2.2 (a) shows a prototype array that is optimized to a given configuration as compared with the source dimensions, Fig. 2.2 (b) shows how to use sequential measurements to make the array denser, and Fig. 2.2 (c) how to make the array larger. Contrary to simultaneous measurements, sequential measurements are faced with the difficulty that phase relationships are lost between consecutive measurements, yet such information is however vital to signal reconstruction in general [START_REF] Oppenheim | The importance of phase in signals[END_REF]. One strategy is to use a set of fixed reference microphones located in the acoustical field which are fully correlated with the sources so as to indirectly measure phase relationships between consecutive positions of the array ; then phases that are lost can be reconstructed by spectral analysis based methods [START_REF] Leclere | Multi-channel spectral analysis of multi-pass acquisition measurements[END_REF] [57] : Conditioned spectral analysis (CSA) and virtual source analysis (VSA) are used, which are based on the average over all passes of a set of reference channels and on the correction of cross-spectral quantities by using transfer functions and averaged reference spectra. A related method that is based on a linear prediction technique is proposed based on the same "fixed references" strategy in Ref [START_REF] Yoon | A method for the efficient construction of acoustic pressure cross-spectral matrices[END_REF]. In Ref [START_REF] Yoon | A method for the efficient construction of acoustic pressure cross-spectral matrices[END_REF], the concept of the rank of a matrix has been employed as the theoretical development. It is proved that the rank of the full cross-spectral matrix equals to the number of uncorrelated acoustic sources and it is used to determine the number of reference microphones.

Discussion : The effectiveness of these approaches is strongly dependent on the availability of high-quality references (high signal-to-noise ratios and high coherence with the source field) in sufficient number to capture the stochastic dimension of the source field (the number of uncorrelated sources it is composed of) ; in particular the number of references should be greater than or equal to the number of uncorrelated sources. Unfortunately, the provision of references incurs an extra cost and precludes some microphones to be used when the user has to face a limited number of tracks in the acquisition system.

Sequential measurements without reference

Kwon and Kim [59] [60] proposed Moving frame acoustic holography (MFAH) to virtually increase the hologram size and the spatial resolution of the hologram by continuously sweeping the sound field with a line array of microphones. They proposed a procedure to handle the Doppler effects on acoustic holograms. Loyau et al. [START_REF] Loyau | Broadband acoustic holography reconstruction from acoustic intensity measurements. i : Principle of the method[END_REF] proposed broadband acoustical holography from intensity measurement (BAHIM) to reconstruct the sources by using near-field acoustical intensity measurement. These quadratic values are measured independently at each point on the hologram and therefore do not require a phase reference signal corresponding to the source. The geometry of the source surface is taken into consideration for an analysis of the phase gradient function on the basis of the elementary orthogonal components associated with this geometry, which is integrated to obtain the pressure field phase. Antoni [START_REF] Antoni | Full-field reconstruction from scanned measurements without references : the latent variable approach[END_REF] proposed a solution to reconstruct the full source field of the radiating object from sequential measurements, yet without the need for fixed references and without the edge effects of the patch approach. Under the stationary assumption, the problem is shown to boil down to the factorization of a structured covariance matrix. This is undertaken within a Bayesian probabilistic approach, where the source field is encoded by unobservable latent variables which are iteratively reconstructed by using an Expected-Maximization (EM) algorithm. The method is shown to return virtually similar results as if all data were captured simultaneously instead of sequentially. In addition, it is ultimately able to simulate statistically equivalent realizations of the source field. A solution is also provided for automatically setting the regularization parameter, an important issue of inverse acoustics. The proposed method does not require any external reference (although it may easily be extended to account for such a possibility) and does not assume the radiated field to be entirely covered by the array. The only assumption it relies on is the stationarity of the source field. This is probably a rather weak assumption, but fundamental for measurements taken at different times to give different images of the same source distribution. Indeed, under the stationary assumption, the problem is shown to boil down to the factorization of a structured covariance matrix. The latent variable model is not only efficient from a computational point of view, but it also enjoys a physical interpretation : the latent variables are "modal coordinates"of underlying virtual incoherent sources that produce the sound field. Discussion : The above EM method directly solves the inverse problem with the incomplete spectral matrix, and this can be seen as an intrusive method, which is a different approach than that followed in this thesis. The methods proposed in this thesis aims to complete a spectral matrix, then feed it to a classic inverse algorithm, which can be understood as non-intrusive. The difference between the two methods are shown in Fig. 2

Conclusion of the chapter

Three categories methods were surveyed in this chapter : Single measurement, Sequential measurements with references, and Sequential measurements without reference. A representative approach for Single measurement is Patch NAH, based on the assumption that the measured signal is band-limited which in turn implies a high correlation of the acoustical field. This is an efficient and easy-to-implement method, however, it suffers from residual prediction errors that increased with the area of the extrapolation surface ; another representative approach for Single measurement is compressed sensing NAH, based on the assumption that the measured signal is sparse under some over-complete dictionary. This method is theoretically guaranteed to be successful with a high probability, however, it is limited to specific surfaces. A representative approach for Sequential measurements with references is based on minimum mean square error linear prediction, based on the assumption that the measured signal is stationary. This is an efficient method but it is strongly dependent on the availability of high-quality references and it also incurs extra costs. Two approaches may be identified for sequential measurements without references. The first one is intrusive, which requires reformulating the inverse problem to account for missing entries in the spectral matrix. The second one is non-intrusive, which is proposed in this thesis. Non intrusive methods aims at directly completing the spectral matrix and then feed it to a classic inverse method.

Setup

Representative

Structured low rank model and Cyclic Projection

In this chapter, a structured low rank model and its corresponding algorithm "Cyclic Projection" are introduced. First, forward problems with sequential and simultaneous measurements are introduced in section 3.1 by discretizing Fredholm integral equations of the first kind. The difference between sequential and simultaneous measurements is then explained as a problem with missing entries in the spectral matrix. In section 3.2, this problem is formulated by means of a structured low rank model : the aim is to find a full spectral matrix subject to given constraint of hermitian symmetry, measurements fitting, reduced rank, and spatial continuity of the sound field ; its ingredients are explained in separate subsections. A Cyclic Projection algorithm with predefined sets and projection operators is proposed to solve the problem in section 3.3.

Stochastic modeling of the forward problem

Stochastic modeling of sequential measurements

Let s(r, ω; ζ), r ∈ Γ be the normal velocity or acoustic pressure on the source surface Γ enclosing the radiating object D, at a given frequency ω. The fundamental premise is that the acoustic source distribution is considered as a stationary stochastic field (only the temporal stationarity is assumed), thus the acoustic sources produce outcomes whose realizations depend on the events ζ in sample space Ω. From a practical point of view, an outcome of the source distribution will simply correspond to a snapshot of the measured signal (i.e. the Fourier transform of a short-time segment possibly tapered with a smooth data window) at a given position of the array.

The acoustic field that is produced by the radiating body D is measured at some discrete locations r m,i , where r m,i denotes the position of m-th microphone in the array, m = 1, ..., M at the i-th position of the array, i = 1, ..., P (an array of M microphones with P measurements), and p(r m,i , ω; ζ) is the acoustic pressure measured at this point [START_REF] Dong | Blind separation of sound sources from the principle of least spatial entropy[END_REF] for technique details to choose the parameters and some applications in related works. Note that for a given array position i, M measurements are taken simultaneously. Similarly, for a given microphone m, the P measurements that result from moving the array at different positions are sequential ; this implies that phase relationships between measurements at points p(r m,i , ω; ζ) and p(r n,j , ω; ζ), ∀(n, m), i ≠ j are necessarily lost.

Having defined the measured pressures and the source field, the direct problem typically relates them through an integral equation of the form

p(r m,i , ω; ζ) = Γ G(r m,i , ω r)s(r, ω; ζ)dΓ(r) + n(r m,i , ω; ζ), ( 3.1)
where G(r m,i , ω r) denotes the Green function (provided in analytical or numerical form) between r and r m,i , and where n(r m,i , ω; ζ) stands for additive measurement noise statistically independent of s(r, ω; ζ). In order to discretize the problem, let us project it onto a basis of spatial functions {φ k (r, ω)} K k=1 whose choice may be either arbitrary or optimal given the object and array geometries :

s(r, ω; ζ) = K k=1 c k (ζ)φ k (r, ω). (3.2)
For a given frequency ω, coefficients c k (ζ) assigned to the basis functions φ k (r, ω) are to be interpreted as K random variables that produce the stochastic field s(r, ω; ζ). The basis of spatial functions can be chosen differently with the various practical setups, and this has been surveyed in previous chapter. One optimal basis is constructed in Ref [START_REF] Antoni | A bayesian approach to sound source reconstruction : Optimal basis, regularization, and focusing[END_REF] : it declares the optimal basis functions are the eigen-functions of a specific continuousdiscrete propagation operator, with the number equals to the number of microphones in the array, which is used in this thesis. Equation (3.1) then becomes

p(r m,i , ω; ζ) = K k=1 H m,i,k c k (ζ) + n(r m,i , ω; ζ), (3.3)
where

H m,i,k = Γ G(r m,i , ω r)φ k (r, ω)dΓ(r)
is interpreted as a transfer function. Equation (3.3) can be recast into a more compact form. let N i be the number of snapshots available at position i , p ij ∈ C M the column vector containing the M pressures {p(r m,i , ω; ζ)} M m=1 measured at the j-th snapshot and i-th position of the array,

c ij ∈ C K the column vector containing the K coeffi- cients {c k (ζ ij )} K k=1 , n ij ∈ C M the column vector of additive noises {n(r m,i , ω; ζ ij )} M m=1 and H i ∈ C M ×K the matrix with elements [H i ] mk = H m,i,k . The matrix version of Eq. (3.3) then reads p ij = H i c ij + n ij , m = 1, ..., M, i = 1, ..., P, (3.4)
where n ij has covariance matrix β i 2 Ω i with assumed spatial coherence Ω i normalized such that its trace tr{Ω i } = M and where the noise power β 2 i is an unknown hyper-parameter. Let us define the expected value of any function f (i) , evaluated at the i-th position of the array, as

E{f (i) } = lim N →∞ 1 N N j=1 f (i) (ζ ij ). (3.5) Next, let S (i) pp ≐ E{p ij p * ij }
(with * the transpose-conjugate operator) define the spectral matrix of the measurements at position i, whose diagonal (resp. off-diagonal) elements contain the auto (resp. cross) spectra of the measured pressures. Similarly, let S cc ≐ E{c ij c * ij } be the covariance matrix of the unknown coefficients, which does not depend on index i by definition of the source field stationarity. Let us also introduce

Ŝ(i) pp = 1 N i N i j=1 p ij p * ij , (3.6)
an estimate of the theoretical spectral matrix S (i) pp obtained by averaging the outer-product of measured pressure vectors over a finite number of snapshots. Therefore, it results from Eq. (3.6) that

Ŝ(i) pp ≃ H i S cc H * i + β i 2 Ω i , i = 1, ..., P, (3.7)
with exact equality when the number of snapshots tends to infinity, N i → ∞, and S cc is a K × K matrix. These covariance equations physically reflect transfer of energy from the source field to the microphones.

Stochastic modeling of simultaneous measurement

In the simultaneous measurement, Eq. (3.4) can be directly written as follows (without subscript i since measurements are captured simultaneously) in a vector form :

p j = Hc j + n j .
(3.8)

The theoretic spectral matrix S pp of simultaneous measurement is then

S pp ≐ E{p j p * j } = HE{c j c * j }H * + E{n j n * j }, where E{n j n * j } = σ 2 n Ω and the definition of S cc ≐ E{c ij c * ij } = E{c j c * j }
remains the same as in the sequential measurements owing to the stationarity assumption of source field. Similar to sequential measurements, Ŝpp = 1 N ∑ N j=1 p j p * j ; then Eq. (3.8) can be written as :

Ŝpp ≃ S L pp + σ 2 n Ω = HS cc H * + σ 2 n Ω, (3.9)
where σ 2 n is the noise power in simultaneous measurements, trace(Ω) = M P , and Ω is a diagonal matrix by assuming that the noise in different channel is uncorrelated. The equality is exact when the number of snapshots N tends to infinity. It is assumed that S L pp = HS cc H * is a low rank matrix, where H ∈ C M P ×K is full row rank and rank(HS cc H * ) = rank(S cc ).

Problem statement

An essential question to be addressed now is the difference between sequential measurements and simultaneous measurement from the spectral matrix perspective. With sequential measurements, an incomplete spectral matrix Ŝm pp is constructed by rearranging all spectral sub-matrices Ŝ(i) pp , i = 1, ..., P in Eq. (3.7) (each with size M × M ) methodically in block diagonal positions and the remaining positions are padded by zero elements (see left graph in Fig. 3.3, it is noted that zero padding is a just a way to visualize the measured matrix and it is in fact a unknown part). Contrary to sequential measurements, Ŝpp in Eq. (3.9) is a complete matrix as shown in the right graph of Fig. 3.3. With the purpose of using sequential measurements for acoustic source reconstruction, the problem is formally stated as to attain an full spectral matrix from the data missing spectral matrix Ŝm pp . In fact, the spectral matrix in sequential measurements with references is also a data missing matrix but with known cross-spectra between the references and each sequential measurements, which is shown as a narrow band in the middle graph of Fig. 3.3. Thus, when the problem is formulated as matrix completion, it is general enough to cover the two cases with or without references.

Structured low rank model

In order to solve the matrix completion problem, the following structured low rank model is proposed : A full spectral matrix S is sought under four constraints which are explained separately hereafter.

find S ∈ C M P ×M P such that ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ Rank(S) = r A (S) -Ŝm pp F ≤ ΨSΨ * -S F ≤ ε S * = S ⪰ 0.
1) Low rank spectral matrix : one wants to find matrix S with given rank r which estimates the low rank S L pp = HS cc H * matrix in Eq. (3.9). This is explained in section 3.2.1.

2) System equation : the system equation A (S) + E = Ŝm

pp represents the data fitting relation, where the unknown spectral matrix S to be searched is linked to the partial measurements Ŝm pp by the sampling operator A ∶ C M P ×M P → C M P ×M P that extracts the elements from diagonal blocks, and where E is composed of estimation and measurement errors. This system equation is imposed as a constraint A (S) -Ŝm pp F ≤ in the model (3.10), where the difference between A (S) and Ŝm pp in Frobenius norm is less than a given tolerance .

3) Spatial continuity of the acoustic field : in order to ensure spatial continuity of the acoustic field, the information on microphone positions must be encoded in the spectral matrix : ΨSΨ * -S F ≤ ε, where Ψ is a projection basis ; this is explained in section 3.2.2.

4) Hermitian property and positive semi-definiteness of the spectral matrix : S * = S ⪰ 0, where ⪰ denotes positive semi-definite. The hermitian property is an inherent characteristic of spectral matrices, as well as the fact that its eigenvalues are nonnegative.

It is noted that the third constraint makes the problem at hand rather different than other matrix completion problems found in the literature [START_REF] Candès | Exact matrix completion via convex optimization[END_REF] [START_REF] Boyd | Alternating projections[END_REF], which is the reason why the proposed model is named as "structured low rank" model. The low rank prior and spatial continuity of the acoustic field are explained in the following subsections. In acoustics, the rank of the spectral matrix has its specific physical implications. The rank of S L pp measures the stochastic model complexity of the acoustic field that is the number of latent variables (uncorrelated sources) required to produce the pressure measurement p [START_REF] Kompella | A technique to determine the number of incoherent sources contributing to the response of a system[END_REF]. The conception of latent variables are introduced : the independent variables with independent and identically distribution (iid), standardized complex random variables l , l = 1, ..., n s , n s < K are used to represent the coefficients c k (ζ) in Eq. (3.3) as

Low rank spectral matrix 3.2.1.1 Low rank in acoustics

c k (ζ) = ns l=1 λ kl l (ζ), k = 1, ..., K, (3.11)
and its matrix form is

c j = Λ j , j = 1, ..., N j , ( 3.12)
where Λ ∈ C K×ns and j ∈ C ns . With this matrix form equation, the Eq. (3.8) is reformulated as

p j = HΛ j + n j , j = 1, ..., N j , ( 3.13)
and the corresponding spectral matrix S L pp can be represented as

S L pp = HS cc H * = HE{c j c * j }H * = HΛE{ j j * }Λ * H * = HΛΛ * H * , ( 3.14)
where H ∈ C M P ×K , and E{ j j * } = I ns . It is obvious that the rank of matrix S L pp and S cc is equal to number of latent variables n s which is physically interpreted as the number of acoustic sources. Furthermore, the rank of the spectral matrix indicates the correlation of the acoustical field. For instance, a rank one spectral matrix corresponds to a fully coherent field, whereas a spectral matrix with full rank corresponds to a diffuse field in general. When the spectral matrix S L pp is low rank, there exists a hidden low dimension structure in the spectral matrix : Fig. 3.4 shows an example of a spectral matrix S L pp under a setup with 3 sources. The full spectral matrix S L pp is shown in Fig. 3.4 (a) with corresponding eigenvalues in Fig. 3.4 (b) ; all the energy is concentrated in the first 3 eigenvalues which indicate the number of sources ; the spectral matrix Ŝm pp from sequential measurements is shown in Fig. 3.4 (c) which is made only of diagonal blocks (unmeasured parts are filled with zeros to show the image) with corresponding eigenvalues in Fig. 3.4 (d), and the energy is nearly scattered over the whole eigenvalues axis.

A Naive solution

The low rank assumption can be seen as an a priori information to be imposed to solve the inverse problem from a Bayesian perspective [START_REF] Salakhutdinov | Bayesian probabilistic matrix factorization using markov chain monte carlo[END_REF] [71] [START_REF] Zhou | Low-rank modeling and its applications in image analysis[END_REF]. Thus, it shrinks the space of solution by adding a constraint. Accordingly, the spectral matrix completion problem can be formulated as fitting a low rank spectral matrix S to the partial measurements Ŝm pp , i.e.

find S ∈ C M P ×M P such that

⎧ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎩ Rank(S) = r A (S) -Ŝm pp F ≤ S * = S ⪰ 0. (3.15)
Problem (3.15) can be solved by Alternating Projection [START_REF] Boyd | Alternating projections[END_REF] or other methods [START_REF] Candès | Exact matrix completion via convex optimization[END_REF] [START_REF] Candes | Matrix completion with noise[END_REF] under some conditions on the number and the distribution of the elements in the matrix, the coherence of the eigenvectors, etc. Unfortunately, it has no unique solution for block diagonal matrices of the type Ŝm pp . An example is given here to illustrate this fact (see [START_REF] Candès | Exact matrix completion via convex optimization[END_REF] for a proof). Assume that the measured matrix Ŝm pp is a 2-by-2 complex hermitian matrix ⎛ ⎝

1 a -1 a 1 ⎞ ⎠
, where a = a e jϕ is an unknown complex element and the diagonal elements are known measurements without error. By definition the objective matrix S is a rank-1 matrix and it is straightforward to find that a e jϕ = a -1 e -jϕ as a result of hermitian symmetry. Therefore, the solution of Eq. (3.15) is non-unique since it can be written {a ∶ a = 1, ∀ϕ}. Thus, another constraint is needed as explained in the next section.

Enforcing spatial continuity of the acoustic field

The spectral matrix S L pp reflects the correlation between the Fourier coefficients of the measurements, however, it does not contain information on the position of the micro-phones. For example, two adjacent columns in the spectral matrix must become equal if the corresponding microphones are infinitely close to each other (nearly overlapping) due to the spatial continuity of the acoustic field ; accordingly, this information must be incorporated in the current problem. In this work, the spatial position of the microphones is encoded in the matrix completion problem by introducing an additional constraint. It is assumed that the measured pressure p with size M P × 1 can be decomposed onto a dimension reduced spatial basis Φ ∈ C M P ×Kp (K p < M P ) with coefficients ϑ ∈ C Kp×1 as follows :

p = Kp i=1 φ i (r)ϑ i = Φϑ, (3.16)
where φ i (r) is the i-th column of matrix Φ and ϑ i is the i-th entry of vector ϑ. In this thesis, the Fourier basis is chosen, meaning that the acoustical field is decomposed as a sum of "acoustical modes" [START_REF] Photiadis | The relationship of singular value decomposition to wave-vector filtering in sound radiation problems[END_REF] (details as how to construct this basis are given in chapter 5). In other words, a spatial structure is encoded. The coefficients ϑ are related to the measurements as :

ϑ = Φ † p, (3.17)
where † denotes the pseudoinverse of a matrix as a result Φ is generally not invertible.

Then the smoothed pressure p can be represented as [START_REF] Pascal | Patch near-field acoustic holography : Regularized extension and statistically optimized methods[END_REF] p = ΦΦ † p = Ψp, (3.18) where Ψ = ΦΦ † is defined as a projection operator. Now, let S ≐ E{pp * } and S ≐ E{pp * }. Therefore,

S = ΨSΨ * . (3.19)
Equation (3.19) may be rewritten as ΨSΨ * -S F ≤ ε where ε is the representation error between ΨSΨ * and S owing to the Φ truncation. It is important that Ψ imposes a specific structure that encodes the information on microphone positions into the matrix S. This can be seen as another constraint to be added to Eq. (3.15). It is noted that the Ψ is an orthogonal projection operator, which can be decomposed as U ⎛ ⎝

I Kp 0 0 0 ⎞ ⎠ U *
, where UU * = I, I Kp is the identity matrix with n elements, and U is the modal matrix in the eigenvalue decomposition (EVD) of Φ. Then Eq. (3.18) can be written as

p = Ψp = Kp i=1 (p * u i ) * u i . (3.20)
This means that only the p * u i i = 1...K p part corresponding to the eigenvalues equal to one are kept by this projection operator, and this is exactly the property of orthogonal projection operator. It is important to emphasize that the role of the orthogonal projection operator Ψ is to smooth the spectral matrix and to ensure the spatial continuity of the acoustical field.

Cyclic Projection

This section now discusses as how to solve problem (3.10) by Cyclic Projection.

From Alternating Projection to Cyclic Projection

Alternating Projection is a method to find a point at the intersection of two closed convex sets by using a sequence of projections from one set to the other [START_REF] Boyd | Alternating projections[END_REF]. The iteration procedure can be defined sequentially as follows. Suppose that C and D are two closed convex sets (every point on the line segment connecting every pair of points within the set) in R N , and define P C and P D the projections on C and D respectively. The algorithm starts with an arbitrary value x k ∈ C and then generates the iteration sequence :

⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ y k = P D (x k ) x k+1 = P C (y k ), k = 0, 1, 2...
where x k ∈ C and y k ∈ D. In words, x k ∈ C is projected onto y k ∈ D first, and second y k ∈ D is projected back onto set C, and so forth. One of the simplest cases is when C and D are closed convex sets and C ⋂ D ≠ ∅ (the intersection is non-empty), in which case the convergence rate of the iterations is proven to be linear [START_REF] Cheney | Proximity maps for convex sets[END_REF]. When C ⋂ D = ∅ (the intersection is empty), alternating projections yields a pair of points in C and D that have minimum distance x ⋆y ⋆ 2 = dist(C, D), where x k converges to x ⋆ in C and y k converges to y ⋆ in D [START_REF] Combettes | The foundations of set theoretic estimation[END_REF]. When C or D is not a convex set, alternating projection has no convergence proof in general (despite local linear convergence when the strong regularity condition holds [START_REF] Lewis | Local linear convergence for alternating and averaged nonconvex projections[END_REF]), yet it is still a popular heuristic to solve the optimization problem. Cyclic projection (CP) is an extension of the basic alternating projection algorithm that generalizes the number of sets to more than 2. Specifically, it can find a point at the intersection of k > 2 sets by projecting alternatively onto C 1 , then C 2 , . . . , and then C k . In order to use CP to solve the current problem in this thesis, the sets and the projection operators should be defined separately, which this is the objective of the next section.

Predefined sets and operators

Sampling operator

The sampling operator and its corresponding set are early defined according to the specific structure of the spectral matrix resulting from the sequential measurements mode. Sequential measurements produce a block diagonal spectral matrix, thus C 1 is defined as a (affine) set in the set of n × n Hermitian matrices S n with specified block entries and the sampling operator P C 1 is defined such that

B = P C 1 (A) with B ij = ⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ [ Ŝm pp ] ij ij ∈ Υ [A] ij ij ∉ Υ, (3.21)
where [] ij denotes the (i, j) element of a matrix and Υ denotes the positions of fixed block diagonal entries. The sampling operator forces the actual measurements to be kept unchanged in each iteration step.

Eigenvalue truncation operator

C 2 = S n
+ is defined as the set of Hermitian positive semidefinite n × n matrices with reduced rank. Note this is a non-convex set [START_REF] Fazel | Rank minimization and applications in system theory[END_REF] and the corresponding projection

B = P C 2 (A) is defined as ⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ A = ∑ n i=1 σ 2 i u i u * i = ∑ r i=1 σ 2 i u i u * i + ∑ n i=r+1 σ 2 i u i u * i B = ∑ r i=1 σ 2 i u i u * i .
First, matrix A is factorized into n eigen-elements, then only the first r largest eigenvalues (σ 2 i where i from 1 to r) and their corresponding eigenvectors (u i where i from 1 to r) are preserved. The eigenvalue truncation operator can be recognized as a hard threshold [START_REF] Kyrillidis | Matrix recipes for hard thresholding methods[END_REF] [79] to fulfill the low rank matrix approximation of the spectral matrix.

Spatial continuity enforcing operator

The last closed convex set is defined as

C 3 = {B ∈ C M P ×M P ∶ B = ΨAΨ * , ∀A = A * }.
The corresponding projection P C 3 (A) denotes an operator that maps the elements A onto set C 3 by B = ΨAΨ * . The introduction of this operator is to endow the spectral matrix with information on spatial positions of the microphones, meanwhile enforcing the spatial continuity of acoustic field.

Cyclic Projection

The solution of problem (3.10) is assumed to be located at the intersection of the three sets C 1 ,C 2 , and C 3 , and is solved by CP. First, the initial value Ŝ{0} is given (for example, a full zeros matrix) ; in each step, the measurements are assigned by projecting Ŝ{k} onto set C 1 to produce S{k} . Subsequently, matrix S{k} is decomposed into eigen-elements, and a low rank approximation is proceeded by truncation : this projects S{k} onto set C 2 to produce Š{k} . Last, spatial continuity is enforced by projecting Š{k} onto set C 3 to produce Ŝ{k+1} . In CP, the order of the cycles is not important, thus the procedure that is described here is only one possibility. For instance, it could also be written as C 3 → C 2 → C 1 and so on. The iteration will be stopped when

[ Ŝm ] ij -[ Š{k} ] ij F [ Ŝm ] ij F
≤ SC, ij ∈ Υ (SC is a constant value that is chosen by user) after projecting onto C 1 and C 2 , then forcing the final solution to be low rank. The proposed CP has only three parameters to tune : the first one is the maximum iteration number M max , the second one is the stopping criterion SC which usually depends on the measurement environments and noise level, and the third one is the rank r which reflects the number of uncorrelated sources in the field. (Enforcing spatial continuity : project from set C 2 onto set C 3 ) 8: go to step 3.

Error analysis of the results

Four kinds of errors are considered : (1) background noise, (2) estimation errors of the spectral matrix, (3) model errors (for example the calibration error of the microphones, the position error of the array etc), and (4) representation errors caused by the smoothing operator in Eq. (3.18). The operator P C 1 and P C 2 bring about the errors (1)(2) and (1) separately, and P C 3 brings about the errors (3)(4). The structured low rank model (3.10) can be reformulated as the projection representation with errors terms (ξ, , ε).

find S ∈ C M P ×M P such that ⎧ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎩ P C 2 (S) -S F ≤ ξ (errors (1)) P C 1 (S) -S F ≤ (errors (1)(2)) P C 3 (S) -S F ≤ ε (errors (3)(4)). (3.22)
This projection representation with error terms may provide more flexibility to apply the CP algorithm in practice. For example, when the errors (1)(2) are too large and errors (3)(4) are acceptable, then the user may choose the stopping criteria according to the constraint P C 3 (S {k} ) -S {k+1} F ≤ ε.

Conclusion of the chapter

The problem of sequential measurements without or with limited number of references has been investigated in this chapter. The main issue is that the phase relationships between consecutive snapshots are missing and result in missing entries in the spectral matrix. The approach proposed in this chapter boils down to a spectral matrix completion problem which is modeled by a structured low rank model and solved by Cyclic Projection. In particular, two of the ingredients : (1) the low rank property of the spectral matrix and (2) the continuity of the acoustic field are highlighted.

• Low rank spectral matrix. The high correlation of the acoustical field implies a low rank of the spectral matrix. The rank of the spectral matrix is consistent with the number of "virtual source", i.e. the number of equivalent uncorrelated sources the acoustical field comprises.

• Continuity of the acoustic field. The high correlation of the acoustical field implies that two columns of the spectral matrix tend to become identical as the corresponding microphone positions become arbitrarily close to each other.

Weakly sparse eigenvalue spectrum model and Fast

Iterative Shrinkage Thresholding Algorithm

In the structured low rank model of the previous chapter, the rank of the spectral matrix was considered to be fixed and known (or at least estimated). This is equivalent to assuming that the number of incoherent sources of the sound field is known. However, this information is rarely available in practice and the estimation of the rank of a matrix is a difficult problem in general. Thus, an alternative model is proposed to tackle this difficulty : the eigenvalue spectrum is assumed to be weakly sparse instead of "exactly sparse", which means that the corresponding spectral matrix is of full rank but with a few dominant eigenvalues. The nuclear norm is then minimized instead of imposing a fixed reduced rank. A Fast Iterative Shrinkage Thresholding Algorithm (FISTA) is proposed to solve the problem. Section 4.1 gives a very short introduction to sparse and low rank constraints before we get into our practical problem. The weakly sparse eigenvalue spectrum model is described in section 4.2. Then, a general mathematical formulation of FISTA is introduced as a preparation to solve the current problem in section 4.3. The FISTA algorithm is finally explained in section 4.4.

A very short introduction to sparse and low rank constraints

This section introduces the basic definitions and conceptions about sparsity and low rank models. With the explosion of massive amounts of high-dimensional data in science and engineering, new tools are required to extract intrinsic structure in high-dimensional data so as to alleviate the curse of dimensionality [START_REF] Markovsky | Low Rank Approximation : Algorithms, Implementation, Applications[END_REF] [START_REF] Candès | Robust principal component analysis ?[END_REF]. Sparse representation is a rapidly evolving field of research [START_REF] Elad | Sparse and redundant representations : from theory to applications in signal and image processing[END_REF], with several applications in acoustics [START_REF] Chu | Robust bayesian super-resolution approach via sparsity enforcing a priori for near-field aeroacoustic source imaging[END_REF] [45] [82] [START_REF] Malioutov | A sparse signal reconstruction perspective for source localization with sensor arrays[END_REF]. A structured signal can be decomposed with a few elements in an over-complete dictionary that reflects the low dimensional structure of the original signal. First, a un-derdetermined system of equations is considered :

y = Ax, (4.1)
where A is a m × n matrix with m ≪ n and x ∈ R n , y ∈ R m . Assume x is k-sparse meaning that the number of nonzeros in the vector x is k, where k = x l 0 ≪ n and x l 0 ∶= #{i ∶ x i ≠ 0} denotes the number of nonzeros in vector x ; then a sparse solution x can be obtained by solving the l 0 norm optimization problem as follows, minimize

x x l 0 subject to y = Ax. (4.2)
To solve problem (4.2) straightly is hopeless since it is a NP-hard problem [START_REF] Bruckstein | From sparse solutions of systems of equations to sparse modeling of signals and images[END_REF]. However it can be relaxed and solved by l 1 norm minimization, which is a smooth approximation (convex envelope) of l 0 . Thus, problem (4.2) can be reformulated as a convex optimization problem [START_REF] Elad | Sparse and redundant representations : from theory to applications in signal and image processing[END_REF],

minimize x x l 1 subject to y = Ax, (4.3)
where x l 1 ∶= ∑ n i=1 x i which is the sum of the absolute values of all components of x. It is noted that if y is contaminated by noise, then Eq. ( 4.3) is reformulated as [START_REF] Ma | Fixed point and bregman iterative methods for matrix rank minimization[END_REF] minimize

x x l 1 subject to Ax -y l 2 ≤ ε, (4.4)
or its Lagrangian version minimize

x λ x l 1 + Ax -y l 2 , ( 4.5)
where l 2 is the Euclidean norm of a vector , x l 2 ∶= ∑ n i=1 x i 2 , ε is a parameter reflecting the level of noise, and λ is a regularization parameter. The model (4.3) can be transformed into a linear programming (LP) problem [START_REF] Antoniou | Practical optimization : algorithms and engineering applications[END_REF] as follows. minimize

x,d n i=1 d i subject to -d i ≤ x i ≤ d i y = Ax, (4.6)
and (4.4) can be represented as a second order cone programming (SOCP) problem [START_REF] Antoniou | Practical optimization : algorithms and engineering applications[END_REF] (in fact, the linear programming problem is a subset of the SOCP problem), which are currently both solved by the interior point method. In order to compute the Newton direction in the interior point method, a linear equation must be solved, which is problematic when the dimension of the coefficient vector and matrix A is large and results in high computational complexity. Thus, one usually resorts to first order methods (i.e. methods that have at most linear local error) to solve the l 1 norm minimization. The Bregman iterative algorithm and iteration shrinkage thresholding (IST) algorithm are two efficient first order iteration methods to solve this problem. Firstly, Bregman iterative algorithm [START_REF] Yin | Bregman iterative algorithms for l1minimization with applications to compressed sensing[END_REF] is inductively defined as follows, starting with x 0 = p 0 = 0.

⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ x k+1 = arg min x∈R n {D p k J (x, x k ) + 1 2 Ax -y l2 } p k+1 = p k -A T (Ax k+1 -y), (4.7)
where the Bregman distance D p k J based on a convex functional J(⋅) between points x and x k is defined as

D p J (x, x k ) = J(x) -J(x k ) -⟨p, x -x k ⟩
, where p = ∇J(x k ) denotes gradient of function and here the function J(⋅) is specifically x l 1 . This Bregman iteration can also be writen in a "adding back the residual" way as follows.

⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ x k+1 = arg min x∈R n {(λ x l 1 + 1 2 Ax -y k l2 } y k+1 = y k + y -Ax k+1 . (4.8)
The above iteration can be further written as the linearized bregman iteration [START_REF] Cai | Linearized bregman iterations for compressed sensing[END_REF] by introducing the idea of fixed point iteration,

⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ x k = shrink(v k-1 , τ ) v k = v k-1 + δ k A T (y -Ax k ), (4.9)
where the shrink operation is defined as shrink(v k , τ ) = max{ v k iτ, 0}sgn(v k i ), δ k is the step size, τ is the shrinkage parameter and sgn is sign function. Another iterative shrinkage thresholding (IST) [START_REF] Elad | Image denoising with shrinkage and redundant representations[END_REF] [89] algorithm can be alternatively to solve the problem (4.4) as follows,

⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ x k = shrink(v k-1 , τ δ k-1 ) v k = x k + δ k A T (y -Ax k ). (4.10)
It is obvious that the IST algorithm is obtained by replacing the v k in the second iteration step of linearized bregman iteration using the x k and the shrinkage parameter are adjusted properly. Further theory about uniqueness of the sparsest solution and performance analysis can be found in Ref. [START_REF] Elad | Sparse and redundant representations : from theory to applications in signal and image processing[END_REF], and other algorithms for finding a sparse solution can be found in [START_REF] Cai | Orthogonal matching pursuit for sparse signal recovery with noise[END_REF] [86] [START_REF] Hale | A fixed-point continuation method for l1regularized minimization with applications to compressed sensing[END_REF].

A "low rank model" can be seen as an extension of sparse representation : the low rank property of a matrix is a natural generalization of the sparsity of a vector [START_REF] Fazel | Matrix rank minimization with applications[END_REF] [START_REF] Recht | Guaranteed minimum-rank solutions of linear matrix equations via nuclear norm minimization[END_REF], which is defined as a matrix with only a few non zero singular values. The low rank model assumes that the information of interest has a compact structure embedded in a high dimensional space [START_REF] Candès | Robust principal component analysis ?[END_REF] [START_REF] Zhou | Low-rank modeling and its applications in image analysis[END_REF], and it has some typical applications in system identification, machine learning, and signal processing [START_REF] Markovsky | Application of structured total least squares for system identification and model reduction[END_REF] [95] [START_REF] Decoste | Collaborative prediction using ensembles of maximum margin matrix factorizations[END_REF] [97] [START_REF] Ji | Robust video denoising using low rank matrix completion[END_REF]. The conceptions of a low rank constraint can be correspondingly represented as follows. First, Eq. (4.1) can be generalized to a matrix equation as follows,

Y = A (X), (4.11)
with the linear map A ∶ R m×n → R m×n and X ∈ R m×n , Y ∈ R m×n . Assume X is a rank-r matrix (where r ≪ min{m, n}), then the low rank solution to the matrix equation (4.11) can be written as minimize

X Rank (X) subject to Y = A (X), (4.12)
Note that Rank(X) is the number of nonzero singular values of X, which is equivalent to the l 0 norm of the singular values vector σ(X)

with σ 1 (X) ≥ σ 2 (X) ≥ ... ≥ σ r (X) > 0 = σ r+1 (X) = ... = σ min{m,n} (
X). Thus (4.12) can be rewritten as minimize

X σ(X) l 0 subject to Y = A (X). (4.13)
Relaxing the l 0 norm to the l 1 norm as before, the problem can be reformulated as minimize

X σ(X) l 1 subject to Y = A (X). (4.14) Let X * ∶= ∑ min{m,n} k=1 σ k (X) = σ(X) l 1
define the nuclear norm, with σ k (X) the kth singular value of X ; then the nuclear norm minimization problem can be written as [START_REF] Cai | A singular value thresholding algorithm for matrix completion[END_REF] [77] minimize

X X * subject to Y = A (X). (4.15)
similarly to (4.4) and (4.5), its version accounted for the presence of additive noise can be separately written as minimize

X X * subject to A (X) -Y F ≤ ε, (4.16)
and in Lagrangian form, minimize

X λ X * + A (X) -Y F , (4.17)
where F is the Frobenius norm of a matrix, ε is a parameter reflecting the level of noise, and λ is a regularization parameter. Further theory can be found in Ref. [67] [100]. Nuclear norm heuristics can be formulated in terms of Semidefinite Programming (SDP) [START_REF] Vandenberghe | Semidefinite programming[END_REF] : minimize

X,W 1 ,W 2 trace(W 1 ) + trace(W 2 ) subject to Y = A (X) ⎛ ⎝ W 1 X X * W 2 ⎞ ⎠ ⪰ 0 (4.18)
where "trace" denotes the trace of a matrix and ⪰ 0 denotes positive semidefiniteness of a matrix, which can be solved by interior point methods [START_REF] Liu | Interior-point method for nuclear norm approximation with application to system identification[END_REF]. The interior point based methods are problematic with the increase of the size of matrix since huge systems of equations need to be solved to calculate the Newton direction. First order methods were proposed to handle this problem. One kind of first order method is the Singular Value Thresholding (SVT) [START_REF] Cai | A singular value thresholding algorithm for matrix completion[END_REF] algorithm, which is inductively defined as follows, starting with

G 0 = 0 ∈ R m×n . ⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ X k = shrink(G k-1 , τ ) G k = G k-1 + δ k (Y -A (X k )), (4.19)
where δ k is the step size, τ is the shrinkage parameter, the shrinkage operation is defined as

shrink(G k-1 , τ ) ∶= Udiag(max{σ i -τ, 0}sgn(σ i ))V * , ( 4.20)
with U, diag(σ i ), V the SVD elements of matrix G k-1 . This iteration can be understood as the linearized Bregman iteration [START_REF] Yin | Bregman iterative algorithms for l1minimization with applications to compressed sensing[END_REF] as explained above. Another kind of first order method, the "iterative shrinkage thresholding" (IST), is proposed in [START_REF] Ma | Fixed point and bregman iterative methods for matrix rank minimization[END_REF], which is inductively defined as follows, (the simulation setup in this example is explained in Fig. 5.1 (b).)

⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ X k = shrink(G k-1 , δ k τ ) G k = X k + δ k (Y -A (X k )), (4.21
In the previous chapter, the proposed structured low rank model was formulated as follows :

find S ∈ C M P ×M P such that ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ Rank(S) = r A (S) -Ŝm pp F ≤ ΨSΨ * -S F ≤ ε S * = S ⪰ 0. (4.22)
The rank of the spectral matrix was considered as a fixed truncation parameter r, which is known a priori or can be estimated with sophisticated statistical method [START_REF] Rissanen | Modeling by shortest data description[END_REF] [START_REF] Akaike | A new look at the statistical model identification[END_REF]. It is noted that the estimation of the rank is an arduous task in general. k (S) is the nuclear norm [START_REF] Fazel | Matrix rank minimization with applications[END_REF] of matrix S computed on the eigenvalues σ 2 1 (S) ≥ σ 2 2 (S) ≥ ... ≥ σ 2 M P (S). In the structured low rank model of the previous chapter, the spectral matrix was assumed to be low rank which implied a sparse eigenvalue spectrum ; a weakly sparse eigenvalue spectrum model assumes the eigenvalue spectrum to be weakly sparse, which involves the l1 norm of the eigenvalues vector of the spectral matrix.

Fast Iterative Shrinkage Thresholding Algorithm

In this section, a fast iterative shrinkage thresholding algorithm (FISTA) [START_REF] Beck | A fast iterative shrinkage-thresholding algorithm for linear inverse problems[END_REF] [108] is introduced to solve Eq. (4.23). First, the general FISTA is introduced.

From projection to proximity operators

The proposed structured low rank model can be generally formulated as a problem to find a point in the intersection of several sets,

find x ∈ m ⋂ i=1 C i where C i is a closed convex set. (4.24)
To find a point x in the intersection of several sets C 1 ... C m ; CP activates each set

C i , i = 1, ..., m
individually by means of its projection operator P C i , and its iteration update rule can be written as where I C (y) is defined as

x k+1 = P C 1 ...P Cm (x k ). ( 4
I C (y) = ⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ 0, if y ∈ C ∞, if y ∉ C
. The projection operator actually finds a point y in C and meanwhile it has minimum distance to x. The sought weakly sparse eigenvalue spectrum will be seen to minimize a combination of objective functions, minimize

x f 1 (x) + f 2 (x) + ...f m (x), (4.27)
with respect to some optimization variable x. A major challenge that arises in solving this problem stems from some of the functions not being differentiable, which excludes conventional smooth optimization techniques. In order to deal with this difficulty, the functions f 1 , ..., f m are considered individually so as to yield an easily implementable algorithm for each nonsmooth function, by proximity operators. Proximity operator prox f is an extension of a projection operator from sets to functions. prox f (x) = y implies finding

y

Fast Iterative Shrinkage Thresholding Algorithm in general

An unconstrained continuously differentiable function minimization problem is considered, and it can be formulated as follows,

minimize x∈R f (x), (4.30)
where function

f ∶ R n → R is smooth convex function with Lipschitz gradient, ∇f (x) -∇f (y) 2 ≤ L f x -y for every x, y ∈ R n , ( 4.31)
where ∇f (x) denotes gradient of function f (x) and L f is Lipschitz gradient constant. It is well known that a smooth convex function f (x) satisfies the following inequality,

f (x) ≥ f (y) + ⟨∇f (y), x -y⟩, (4.32)
where the inner product ⟨∇f (y), x-y⟩ is calculated as ∇f (y) T (x-y). This inequality implies geometrically that the tangent lines of the function denoted by f (y) + ⟨∇f (y), x -y⟩ are underneath the function curve denoted by f (x). Since the function f (x) is assumed with Lipschitz gradient, f (x) is also bounded by a family of local convex quadratic functions with Hessain matrix L f I and thus,

f (y) + ⟨∇f (y), x -y⟩ ≤ f (x) ≤ f (y) + ⟨∇f (y), x -y⟩ + L f 2 x -y 2 2 , ( 4.33)
where f (y) + ⟨∇f (y), x -y⟩ is the linear part of function f at point y and the quadratic proximal term

L f 2
x-y 2 2 measures the local approximation error. Let consider an iteration sequence x k generated by a basic gradient descent algorithm,

x k = x k-1 -µ∇f (x k-1 ), (4.34)
where µ is a step size. According to Eq. (4.33), this gradient iteration can be further viewed as a proximal regularization of the linearized function f at x k-1 as follows, where g(x) is a continuous convex function which is possibly non-smooth and λ is a regularization parameter. Let the scalar step size be µ = 1 L f which is constant in each iteration ; then the solution of Eq. (4.36) can be written iteratively as follows according to Eq. (4.35).

x k = arg min x {f (x k-1 ) + ⟨∇f (x k-1 ), x -x k-1 ⟩ + L f 2 x -x k-1 2 2 }. ( 4
x k = arg min x {f (x k-1 ) + ⟨∇f (x k-1 ), x -x k-1 ⟩ + 1 2µ x -x k-1 2 2 + λg(x)}. (4.37)
Let's complete the square and neglect constant terms ; Eq. (4.37) can then be converted to a square form

x k = arg min x { 1 2µ x -(x k-1 -µ∇f (x k-1 )) 2 2 + λg(x)}. (4.38)
The solution of Eq. (4.38) can be written as a proximal-gradient algorithm (as contrasted to classic gradient algorithm),

x k = prox λµ (x k-1 -µ∇f (x k-1 )), (4.39)
where prox λµ stands for the proximity operator with parameters (λ, µ) and the definition of proximity operator depends on the function g(x). A list of proximity operator for its corresponding function can be found in Ref. [START_REF] Bauschke | Fixed-point algorithms for inverse problems in science and engineering[END_REF]. It is noted that Eq. (4.39) reduces to the gradient descent

x k = x k-1 -µ∇f (x k-1 ) (4.40)
when g(x) = 0 and Eq. (4.39) reduces to proximal point algorithm when f (x) = 0,

x k = prox λµ (x k-1 ). ( 4
F (x) = f (x) + λg(x) subject to x ∈ C, (4.42)
where C is assumed to be a nonempty closed convex set ; then the iterative solution can be correspondingly written as

x k = P C (prox λµ (x k-1 -µ∇f (x k-1 ))), (4.43)
where P C is the projection onto set C. Equation (4.43) shares a sublinear global rate of convergence (the convergence speed of sequence {x k } to the optimal solution x ⋆ )

F (x k ) -F (x ⋆ ) ≃ O(1 k). ( 4.44)
A fast algorithm to accelerate the convergence was first investigated by Nesterov [START_REF] Nesterov | A method of solving a convex programming problem with convergence rate o (1/k2)[END_REF] and later rediscovered by Beck and Teboulle [111] [107]. Therefore, the fast Iterative algorithm can be formulated as follows,

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ y k = y k-1 -µ∇f (y k-1 ) xk = prox λµ (y k ) ≜ shrink(y k , λ k µ) x k = P C (x k ) t k+1 = 1 2 (1 + 1 + 4t 2 k ) y k+1 = x k + t k -1 t k+1 (x k -x k-1 ), (4.45)
where t k is introduced as an instrumental parameter to accelerate the iteration. When the proximity operator prox λµ is calculated as a shrinkage operator prox λµ (y k ) ≜ shrink(y k , λ k µ) which is one kind of proximity operator corresponding to a category of functions g(x), the algorithm is named Fast Iterative Shrinkage Thresholding Algorithm (FISTA) and its convergence rate is shown to be F is a convex function with Lipschiz gradient L f and this term can be recognized as f (x) in the generic model. S * is a non-smooth convex function that can be recognized as g(x) in the generic model. Thus the solution can be written as follows according to Eq. (4.39),

F (x k ) -F (x ⋆ ) ≃ O(1 k 2 ). ( 4
S k = shrink(S k-1 -µ∇f (S k-1 ), λµ), ( 4.49)
where µ = 1 L f and the gradient [START_REF] Cai | A singular value thresholding algorithm for matrix completion[END_REF] (to be more accurate, subgradient [START_REF] Boyd | Convex Optimization[END_REF] here) of the function f (S) = A (S) -Ŝm Let G k = S k-1 -µ∇f (S k-1 ) ; the soft thresholding operater shrink(G k , λ k µ) can be defined as the following optimization problem,

shrink(G k , λ k µ) = arg min Y 1 2 Y -G k 2 F + λ k µ Y * . (4.51)
The soft thresholding operater tries to find a matrix Y as close as possible to G k with minimum nuclear norm, which is calculated as Udiag(max{σ 2 iλ k µ, 0})U * , where U nad diag(σ 2 i ) are the elements of the Eigenvalue Decomposition (EVD) of matrix G k , and max{σ 2 iλ k µ, 0} takes the maximum value between 0 and σ 2 iλ k µ. Thus, according to Eq. (4.45) and with the definition of Eqs. (4.50) and (4.51), the FISTA iteration can be finally formulated as follows

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ G k = G k-1 -µ(A (G k-1 ) -Ŝm pp ) Sk = shrink(G k , λ k µ) S k = Ψ Sk Ψ * t k+1 = 1 2 (1 + 1 + 4t 2 k ) G k+1 = S k + t k -1 t k+1 (S k -S k-1
).

(4.52)

Continuation technique [START_REF] Yin | Bregman iterative algorithms for l1minimization with applications to compressed sensing[END_REF] is used here to tune the regularization parameter λ [START_REF] Hale | A fixed-point continuation method for l1regularized minimization with applications to compressed sensing[END_REF] [84] [START_REF] Toh | An accelerated proximal gradient algorithm for nuclear norm regularized linear least squares problems[END_REF]. It can be seen as a method to solve a sequence of problems (4.36) defined by an increasing sequence {λ k } : λ 0 as an initial regularization parameter and λ d as final regularization parameter are chosen firstly ; a high starting value of λ implies stronger regularization and the optimal manifold [START_REF] Sra | Optimization for machine learning[END_REF] has lower dimension and is easier identify ; by decreasing the regularization parameter λ step by step (with a ratio η) ; a coarse to fine search is achieved in the solution space. It is noted that only a heuristic understanding is given here since there is not much analysis of this approach in the literature. The full algorithm is illustrated in Fig. 4.2 and the whole procedure of FISTA is listed in Algorithm 2. The algorithm is stopped by the following criterion :

Ŝm pp {ij} -G k+1 {ij} F Ŝm pp {ij} F ≤ SC, ij ∈ Υ, (4.53)
where SC is a constant value that is chosen by user.

Conclusion of the chapter

The problem of sequential measurements without reference suffers from the issue that the phase relationships between consecutive snapshots are missing and result in missing entries of the spectral matrix. A Weakly sparse eigenvalue spectrum model has been proposed to solve this problem, which is shown to boil down to minimizing the nuclear norm of a spectral matrix subject to measurements fitting, hermitian symmetry, and spatial continuity of the sound field. A Fast iterative shrinkage thresholding algorithm (FISTA) has been proposed to solve this nuclear norm minimization problem. The method benefits from "weakly sparse" eigenvalue spectrum assumption, which has no requirement of a "fixed rank" as in previous structured low rank model.

Algorithm 2:

Nuclear Norm minimization by FISTA (the Matlab codes are attached in the appendix section)

1: Starts with G 0 = S 0 = 0 ∈ C M P ×M P , t 1 = 1 ; µ is step size, λ 0 is a initial regularization
parameter and λ d is the final regularization parameter, and η is the radio to be decreased for λ k-1 in each step.

2: While λ k ≥ λ d do 3: For 1 ∶ N m 4: G k = G k-1 -µ(A (G k-1 ) -Ŝm pp ) 5: Sk = shrink(G k , λ k µ) 6: S k = Ψ Sk Ψ * 7: t k+1 = 1 2 (1 + 1 + 4t 2 k ) 8: G k+1 = S k + t k -1 t k+1 (S k -S k-

Simulation and experimental validation

The structured low rank model and its corresponding algorithm, Cyclic Projection (CP), as well as the weakly sparse eigenvalue spectrum model and its corresponding algorithm, Fast Iterative shrinkage thresholding algorithm (FISTA), have been explained in Chapter 2 and Chapter 3, respectively. The objective of this chapter is to test CP and FISTA in diverse simulation and experimental setups. A common simulation platform with five kinds of setups is constructed in section 5.1, simulation results about sequential measurements without reference are given in section 5.2, sequential measurements performance analysis with various SNRs is investigated in section 5.2.2, simulation results about sequential measurements with different shift distances are illustrated in section 5.3, results of sequential measurements with references are given in section 5.4, and finally, CP and FISTA are validated with experimental data in section 5.5.

Simulation setups and parameters

In this section, several simulations are constructed to test the sequential measurements algorithms with and without references under various scenarios. Three point sources with equivalent unit magnitudes located separately at (0.2701, -0.0084), (-0.1613, 0.2348), (0.0641, 0.1573) (meters) are simulated to generate the acoustic field. In order to compare the errors of acoustical sources reconstruction hereafter, the source plane is discretized uniformely by a 41 × 41 grid with distance 0.025 m, length 1 m and width 1 m, and these three point sources are expanded in a B-spline basis [START_REF] Antoni | A bayesian approach to sound source reconstruction : Optimal basis, regularization, and focusing[END_REF] to produce a smoothed velocity distribution ; an example of simulated sources is given in Fig. 5.2. The air mass density is 1.2 kg m 3 and the sound velocity is 341 m s. The measurement plane is located at 0.1 m away from the source plane : the microphones are uniformly distributed with spacing 0.1 m in a square array with side length 0.4 m ; the total number of microphones in the array is 25. Complex Gaussian noise is added to 100 snapshots to produce the measured pressures (in frequency domain), and simultaneous measurement spectral matrix (full matrix) and sequential measurements spectral matrix Ŝm pp are separately simulated. All the simulations in this thesis are based upon five setups.

• Setup 1 : simulation of sequential measurements without reference (see Fig. 5.1 (b)). The prototype array is shifted 9 times sequentially at positions (-0.16, -0.16), (-0.16, 0.08), (-0.16, 0.16), (0.08, -0.16), (0.08, 0.08), (0.08, 0.16), (0.16, -0.16), (0.16, 0.08), (0.16, 0.16) (meters). This setup is mainly used to verify and test the CP/FISTA algorithm under diverse choices of bases and SNRs in the case without reference.

• Setup 2 : simulation of sequential measurements with references. The sequential measurements are carried out as in Setup 1, but different numbers of references are considered :

-4 references (see Fig. A performance comparison with the classic method [START_REF] Yoon | A method for the efficient construction of acoustic pressure cross-spectral matrices[END_REF] is provided by using this setup.

Setup 3, 4 and 5 are used to investigate the performance of algorithms with respect to different relative distances between consecutive measurements. The positions of references are as in Setup 2.

• Setup 3 : simulation of sequential measurements without and with references (see Fig. 5.1 (d)) : the prototype array is shifted 9 times sequentially at positions (-0.42, -0.42), (-0.42, 0.08), (-0.42, 0.43), (0.08, -0.42), (0.08, 0.08), (0.08, 0.43), (0.43, -0.42), (0.43, 0.08), (0.43, 0.43) (meters).

• Setup 4 : simulation of sequential measurements without and with references (see Fig. 5.1 (e)) : the prototype array is shifted 9 times sequentially at positions (-0.3, -0.3), (-0.3, 0.08), (-0.3, 0.3), (0.08, -0.3), (0.08, 0.08), (0.08, 0.3), (0.3, -0.3), (0.3, 0.08), (0.3, 0.3) (meters).

• Setup 5 : simulation of sequential measurements without and with references (see Fig. is the completed spectral matrix and (2) the acoustic source reconstruction error (ARE) calculated as

ARE = S {sim} -S {rec} F S {sim} F , (5.2)
where S {sim} is the simulated source and S {rec} is the reconstructed one using the method of Ref [START_REF] Antoni | A bayesian approach to sound source reconstruction : Optimal basis, regularization, and focusing[END_REF] based on the completed spectral matrix.

Sequential measurements without reference

Choice of spatial basis

The purpose of this subsection is to give the guideline to construct the basis, and only CP is used here to verify the results since the basis that are used in CP and FISTA is same. First, the spectral matrix completion results are investigated with the Fourier basis in Setup 1 and the SNR is fixed to 60 dB. The basis functions of the Fourier basis are Φ(x, y) = e i(kxx+kyy) where x, y are the coordinates of the microphones in the array and k x and k y the wavenumbers along the x and y directions respectively (note that Fourier basis is non-regular sampled). The wavenumbers are discretized as k n x = n∆k x (resp. k n y = n∆k y ) n = -N, ..., N where the maximum spatial frequency k max = N ∆k x = π ∆x (resp.N ∆k y = π ∆y ) is related to the minimum spatial resolution ∆x (resp. ∆y) and ∆k x = 2π Lx (resp. ∆k y = 2π Ly ) is the spatial resolution related to the aperture L x (resp. L y ). The first group of bases is constructed with L x = L y = 1 m and ∆x (= ∆y) is set to 0.12, 0.11, 0.10, 0.09, 0.08 (meters), resulting in dimensions of 225×64, 225×81, 225×100, 225×121, 225×144 respectively. The corresponding MCE curves up to 3 kHz are shown in Fig. 5.3 (a). When more and more basis components are introduced up to a certain point about 144 components corresponding to ∆x = 0.08 m, the MCE decreases. However, this increase can not be unlimited (the guideline and the safe bound to construct the basis will be explained later in this section) . The second group of spatial bases is constructed with a measurement aperture L x = L y = 2 m and ∆x is set to 0.24, 0.22, 0.20, 0.18, 0. It is noted that one fundamental assumption in this thesis is that the acoustical field is highly correlated, this assumption implies the correlation length l c of the acoustic field must be larger than the average distance d c between the microphones of sequential measurements : l c > d c (one example is that two sequential measurements can not be carried out too far from each other). Thus, the basic guideline to construct the Fourier basis is as follows : ∆x(∆y) should be larger than the correlation length l c of the acoustic field. Besides, the Nyquist-Shannon sampling theorem requires that ∆x(∆y) should be at least twice as large as the average distance d c between the microphones of sequential measurements. Thus, the rule ∆x(∆y) ≥ min{l c , 2d c } may be used as illustrated in Fig. 5.4. In the following, the basis is constructed with settings ∆x = 0.08 m and L x = L y = 1 m (the corresponding dimension is 225 × 144), although this choice may not be optimal for all the setups.

Discussion : another way to construct the basis is by introducing the analytical Green function, yet this requires certain information about the a priori source distribution ;

Therefore the Fourier basis makes the methods proposed in this thesis truly flexible in the context of a non-intrusive approach. Setup 1 is used here to verify the performance of spectral matrix completion with respect to various values of SNR. Figure 5.5 (a) illustrates the spectral matrix completion results by CP. It is seen that the MCE decreases with the SNR ; although the performance of CP is influenced by noise, the MCE in the whole range is less than 0.12. The phenomenon that MCE in 30 dB are lower than in 60 dB at some frequency range is observed in the simulation is caused mainly by the non optimal choice of the basis function : when the SNR is high (for example SNR = 60 dB and 30 dB), the model errors may play a key role in the final results (the basis here may be a better basis for 30 dB than 60 dB, thus 30 dB case has less model errors) ; when the SNR is relatively low (for example SNR = 20 dB and 10 dB), the noise may play a key role in the final results. The completed spectral matrix is further used to reconstruct the sources by back-propagation [START_REF] Antoni | A bayesian approach to sound source reconstruction : Optimal basis, regularization, and focusing[END_REF]. Figure 5.5 (b) illustrates the acoustic source reconstruction results. It is remarked that the ARE decreases with frequency in general. The spatial aliasing frequency of the prototype array is 1700 Hz, and the ARE by using sequential measurements is far less than 0.01 even at 3000 Hz, which greatly improves the maximum working frequency. It is also observed that the ARE can be quite low even with a relatively high MCE of the spectral matrix ; take for example 3000 Hz with SNR = 10 dB : the MCE is 0.1081, whereas the ARE is 0.0399. This can be explained as a result of regularization which acts as a low pass filter that filtrates errors existing at high frequencies.

Sequential measurements performance analysis with various SNRs

Similar results by FISTA are given in Fig. 5.6, which is competitive with CP. Note that the ARE returned by FISTA in Fig. 5.6 (b) is slightly different from the ARE returned by CP in Fig. 5.5 (b) since matrix completion and source reconstruction are two independent inverse problems as mentioned before ; when the regularization parameters used in source reconstruction based on the spectral matrix completed by CP and FISTA are set identically, then similar results are obtained in Fig. 5.5 (b) and 5.6 (b).

Sequential measurements with different shift distances

The simulation of this section are based on Setup 1 to 5 without reference. Figure 5.7 (a) and Fig. 5.7 (b) show the MCE and ARE respectively, with the shift distances increased from Setup 5 to Setup 4 to Setup 3. It is seen that the MCE decreases with the shift distance. Small shift can capture more correlation of the acoustical field ; on the contrary, if sequential measurements are carried out far from each other, spatial correlation is lost which renders the spectral matrix harder to be completed. Thus, a good spectral matrix completion and acoustical source reconstruction performance depends on the shifts between sequential measurements. problem in setup 3 with the current parameters, where the SVD does not seem to converge after certain steps (this can however be fixed by modifying the parameters). compared with the classical reference-based which is optimal in the mean square sense (Mean Square, MS) method [START_REF] Yoon | A method for the efficient construction of acoustic pressure cross-spectral matrices[END_REF]. The following cases are identified :

Sequential measurements with references

1) The number of references is larger than the number of sources (i.e 4 references).

When SNR = 60 dB, the spectral matrix completion error by MS is far less than by CP/FISTA (see Fig. . When the number of references is greater than or equal to the number of sources in the absence of noise and the MS method is optimal (in the mean square sense) ; the reference signals then provide enough information on correlation between sequential measurements. However the performance of MS rapidly deteriorates when the level of noise increases.

2) The number of references is equal to the number of sources (i.e 3 references). When SNR = 60 dB, the spectral matrix completion error by MS is still far less than by CP/FISTA (see Fig. 3) The number of references is less than the number of sources (i.e 2 references or 1 reference). The MS generally fails to deal with this case, whereas CP and FISTA still work very well. The MCE is illustrated in In conclusion, it is seen that CP/FISTA returns reasonably low errors in all investigated scenarios, contrary to MS which rapidly fails when the number of references is insufficient or becomes noisy. It is remarkable that CP/FISTA achieves more or less the same figure of merit independently of the number of references. Since the CP/FISTA is designed to solve a matrix completion problem, and the matrix completion performance mainly depends on the number of known elements (the more the better) and the rank of the matrix (the small the better) ; and when the rank is small enough (rank equal to three as compared the dimension of matrix is 225 in the simulation), each reference only increase 2 lines known elements about the matrix (as in Fig. 3.3), this insignificantly increased known elements may be difficult to enhance the matrix completion performance significantly. Lastly, a mean MCE and ARE of all the frequencies and SNRs (it calculates the mean of all the curve in Fig. 5.11, 5.12, 5.13 and Fig. 5.14 along the frequency and SNRs, which returns one element for each of number of references) with respect to the increase of number of references are shown in Fig. 5.15, providing a global perspective to the performance of CP/FISTA and MS (note that there exists huge mean error difference between the CP/FISTA and MS since MS is too sensitive to the noise especially when SNR = 10dB).

As an endnote to this section, it should be highlighted that spectral matrix completion and acoustical source reconstruction reflect two different problems. ARE in Fig.

Experimental validation

This section now illustrates the methodology on measurements made in a semi-anechoic chamber. The experimental setup is shown in Fig. 5.16 (a). The sources are 3 Fostex 6301B speakers with 0.1 m diameter and working frequency is from 80 Hz to 13 kHz, plates with absorbent material are placed on the ground in order to limit the reflections of the ground. A rectangular plane array is used with 5 × 6 GRAS microphones P Q40 regularly spaced every 0.1 m. A data acquisition system OROS with 32 channels associated with the software NVgate is used. The proposed experiment makes use of 3 independent sources. Two are generated with the software NVgate (OROS), and the third one is generated externally. The three sources are placed in a plane parallel to the array. The height of the center of the array is 1.17 m. The centers of the two bottom sources are aligned with the center of the array in the initial setup. The center point of the array is considered as the origin of coordinates (0, 0) m. The distance between the source plane and array is 0.73 m and the array is moved at (0, 0), (-0.08, 0), (0.08, 0), (0.08, -0.08), (0, -0.08), (0, 0), (-0.08, -0.08), (-0.08, 0.08), (0, -0.08), (0.08, 0.08) (meters). The sampling frequency is F s = 25.6 kHz. The signals generated for the three sources are independent white noises with frequency bandwidth between 50 Hz and 10 kHz. Two reference signals are acquired directly by two fixed microphones and the third one is taken directly from one signal generator. This provides ideal references with a very high signal to noise ratio. A data file that contains 32 signals (29 microphones and 3 references) are recorded, each with an acquisition time of 10 seconds.

The parameters of CP and FISTA are chosen as same as in the simulation part. The Fourier basis is discretized with ∆x = ∆y = 0.1 m over a surface with length Lx = 1 m and width Ly = 1 m (the dimension of the basis Φ is 225 × 100). First, the CP/FISTA and MS methods are used to complete the spectral matrix in the frequency range from 3000 Hz to 6000 Hz with steps of 100 Hz ; CP is applied without reference and MS is applied with 3 references. The MS results are as a point of comparison to validate the CP/FISTA method, yet with no guarantee that they will return the true values, as shown in the simulation section. The relative completion errors between CP/FISTA and MS methods are shown in Fig. 5.16 (b). It is seen that MCE increases with frequency, while the relative error of acoustical sources reconstruction is lower than 0.3151 in the whole frequency range. The results show that CP/FISTA without reference can achieve comparable performance as the MS (with three references) for acoustical source reconstruction even though the matrix completion results may differ. Figures 5.17 5.18, 5.19 and 5.20 show the acoustical source reconstruction results at 3000 Hz, 4000 Hz, 5000 Hz, 6000 Hz, respectively. It is seen that for a fixed position of the array, the reconstructed sources suffer from limited spatial resolution at low frequencies and their magnitudes are seriously underestimated at high frequencies. Incidentally, the allowable maximum frequency resolution is 1700 Hz. The (k) in Fig. 5.17 -5.20 shows the source reconstruction result by CP. The (l) in Fig. 5.17 -5.20 shows the source reconstruction result by FISTA ; it is seen that at least two sources are well reconstructed even at a high frequency 6000 Hz which greatly extend the maximum working frequency for a fixed position of the array by CP and FISTA. The results are comparable with the MS method. Figure 5.21 shows the spectral matrix completion results at different frequencies with the different methods.

Conclusion of the chapter

CP and FISTA have been validated both by simulation and experimental data, the main results obtained are highlighted hereafter :

1) CP and FISTA can be both used to realize sequential measurements without and with references.

2) CP and FISTA can be considered as a better alternative to MS in general.

Exploring the differences between CP and FISTA

In this thesis, CP and FISTA are investigated for the common goal of completing the spectral matrix ; this can be understood in an uniform way from a signal reconstruction perspective as described in section 6.1. The differences between CP and FISTA are then explored in section 6.2, and eventually simulation and parameter analysis are given in section 6.3.

Understanding CP/FISTA from a signal reconstruction perspective

The upper graph in Fig. 6.1 shows the case where the spatial signal is measured simultaneously without losing the phase information, with the black dots indicating the positions of the microphones. The middle graph shows two sequential measurements, with the black dots indicating the first measurement and the red dots the second measurement ; it is observed that the phase information is lost between the two consecutive measurements. Nuclear norm minimization (or fixed low rank) plays a role to correct the phase of the two consecutive measurements. However, phase connection is still not enough to guarantee uniqueness, the reason why the signal is expanded onto a predefined basis of functions as shown by the blue curve in the bottom-right of Fig. 6.1 : information on the array position is thus encoded and the signal is further corrected and smoothed by the projection operator ΨSΨ * ≈ S.

Comparison between CP and FISTA

Let S be the spectral matrix to recoverer ; it can be decomposed into two parts as follows

S = S lowrank + σ 2 n I, (6.1)
where I is the identity matrix and rank(S lowrank ) = r. Then the EVD (Eigenvalue Decomposition) of the S can be written as

S = UΛU * + σ 2 n I = U(Λ + σ 2 n I)U * = UΠU * , ( 6.2)
where

Π = Λ + σ 2 n I = diag(σ 2 1 + σ 2 n , ⋯, σ 2 r + σ 2 n , σ 2 n , ⋯, σ 2 n ) and Λ = diag(σ 2 1 , ⋯, σ 2 r , 0, ⋯, 0), σ 2 1 ≥ σ 2 2 ≥ ⋯σ 2 r , U = [u 1 , u 2 , .
.., u M P ] is an unitary matrix. The spectra of eigenvalues are illustrated in Fig. 6.2. If the SNR is assumed high enough so that σ 2 r ≥ σ 2 n , then the first r principal eigenvalues

λ 1 = σ 2 1 + σ 2 n , λ 2 = σ 2 2 + σ 2 n , ⋯, λ r = σ 2 r + σ 2 n dominate the remaining M P -r eigenvalues λ r+1 = σ 2 n , λ r+2 = σ 2 n , ⋯, λ M P = σ 2 n
; therefore spectral matrix S can be rewritten as

S = (U s , U n ) ⎛ ⎝ Λ s 0 0 Λ n ⎞ ⎠ ⎛ ⎝ U * s U * n ⎞ ⎠ = U s Λ s U * s + U n Λ n U * n ( 6.3)
where Λ s = diag(λ 1 , ⋯, λ r ) and Λ n = diag(λ r+1 , ⋯, λ M P ), U s = [u 1 , ⋯, u r ] spans the signal subspace and U n = [u r+1 , ⋯, u M P ] spans the noise subspace. It is found that the rank of the spectral matrix is "approximately low", which means the first r eigenvalues are prominent. CP returns an estimate Ŝcp of the signal subspace part,

Ŝcp ≈ U s Λ s U * s , ( 6.4)
which is illustrated at the lower left corner of Fig. 6.2. FISTA assumes that the eigenvalues of matrix S are weakly sparse, which returns an estimate Ŝfista with minimum sum of the eigenvalues,

Ŝfista ≈ U s Λ s U * s -S r + S b , ( 6.5)
where S r is subtracted by virtue of the part that suppresses the noise in U s Λ s U * s , and S b represents the noise part that is introduced by U n Λ n U * n . This is illustrated in the lower right corner of Fig. 6.2. The performance of CP and FISTA can be considered by distinguishing the following cases.

• S r F = S b F indicates that the suppression of noise is equal to the introduction of noise : FISTA will achieve the comparable performance with CP.

• S b F > S r F indicates that the introduction of noise is larger than the suppression of noise : CP will outperform FISTA in this case.

• S b F < S r F indicates that the suppression of noise is larger than the introduction of noise : FISTA will outperform CP in this case.

Although the above analysis provides an intuition to understand CP and FISTA, since S r and S b are influenced by the parameters of the algorithms, and it is hard to conclude whether CP or FISTA is better in general. However, a guideline will be discussed in the next section about as how to choose between CP and FISTA in practice based on the heuristic analysis of this section. Table 6.1 provides a list of the characteristics of CP and FISTA which may help the reader to choose an algorithm in practice. In a nutshell, CP is based on three simple iterations and is quite easy to implement. However it has no guarantee to find the global minimum and the number of acoustic sources needs to be known a priori or estimated. The converge rate is linear [START_REF] Cheney | Proximity maps for convex sets[END_REF] (see more discussions in section 3.3). On the other hand, FISTA is a more sophisticated algorithm that converges faster, like O(1 k 2 ), to the global minimum (convex cost function), yet it has many parameters to tune ; it is not necessary to know the number of sources. Finally, the rank truncation operator (for fixed rank in the structured low rank model) and shrink operator (for minimization of the nuclear norm) can be understood as eigenvalues filter 1 and 2 which is illustrated in Fig. 6.3. Rank truncation operator truncates the eigenvalues spectrum by only keeping the principal eigenvalues (a truncation with hard cut off), while shrink operator can be seen as an eigenvalues filter with smooth transition zone.

characteristics

CP and FISTA with different numbers of sources and

SNRs

The objective of this subsection is to investigate CP and FISTA results with different numbers of sources and SNRs. Setup 1 in chapter 4 is used here to simulate sequential measurements with 3, 9, 16, 25, and 36 sources separately. The parameter settings of the algorithm remain unchanged as in chapter 4 for CP and FISTA, and the analysis of the results is given hereafter.

• Number of sources is 3 : when SNR = 60 dB, CP outperforms FISTA in the frequency range above 1600 Hz and FISTA outperforms CP in the frequency range below 1600 Hz. CP outperforms FISTA when SNR = 30 dB, SNR = 20 dB, SNR = 10 dB as seen in Fig. 6.6 (a), which implies that the introduction of noise is larger than the suppression of noise for FISTA in this situation. Thus, CP may be a better choice when the number of sources is small (CP enforces the rank to be equal to the number of sources) and this number can be estimated by trial and error cheaply, even though the performance of FISTA is also absolutely acceptable in this situation.

• Number of sources is 9 : FISTA outperforms CP for SNR = 60 dB and SNR = 30 dB as seen in Fig. 6.6 (b), which implies that the introduction of noise is less than the suppression of noise for FISTA in this situation. When SNR = 20 dB, CP outperforms FISTA in the frequency range below 1600 Hz and FISTA outperforms CP in the frequency range above 1600 Hz. Similar results are obtained with SNR = 10 dB : CP outperforms FISTA in the frequency range below 2000 Hz and FISTA outperforms CP in the frequency range that is greater than 2000 Hz.

• Number of sources is 16 : FISTA outperforms CP when SNR = 60 dB and SNR = 30 dB as seen in Fig. 6.6 (c) ; when SNR = 20 dB, CP outperforms FISTA in the frequency range below 700 Hz and FISTA outperforms CP in the frequency range above 700 Hz ; when SNR = 10 dB, CP outperforms FISTA in all the frequency band which it implies that the introduction of noise is larger than the suppression of noise for FISTA in this situation. This suggests using CP instead of FISTA when the SNR is small. It will also be found in the next section that the CP performance is not influenced too much by the underestimation and overestimation of sources.

• To sum up, some empirical guidelines to choose the algorithm are listed in Table 6.2. When the SNR is high, it is ealier to have the introduction of noise less than the suppression of noise for FISTA, and vice versa for CP when the SNR is low. When the number of sources is small, CP seems to outperform FISTA. For large numbers of sources, FISTA seems preferable. criteria algorithm high SNRs FISTA low SNRs CP number of sources is large FISTA number of sources is small CP Table 6.2 -Guideline to choose the algorithms

Parameter analysis for CP

The above analysis have been carried out by assuming that the number of sources is correctly known/estimated, which implies that best results of CP are compared with FISTA. In the subsequent subsection, the effects of overestimation and underestimation of the number of sources in CP are investigated.

Overestimation of the number of sources

Figure 6.7 investigates the overestimation by 3 sources (denoted by CP+3) and compares it to correct case (denoted by CP). When the number of sources is 9, CP is better than CP+3. When the number of sources is increased to 16, 25 and 36, the results of CP is comparable with CP+3 since the 3 extra sources that are accounted for corresponds to small eigenvalues. Figure 6.8 compares the overestimation by 3 sources (denoted by CP+3) and 6 sources (denoted by CP+6). It is found that all the results are comparable, which means that CP is robust to the overestimation of the number of sources.

Underestimation of the number of sources

Figure 6.9 compares the underestimation by 3 sources (denoted by CP-3) and compares it to the correct case (denoted by CP). When the number of sources is 9, CP-3 is better than CP, which implies that the introduction of noise is less than the suppression of noise for CP. When the number of sources increased to 16, 25 and 36, the results of CP is comparable with CP-3 since the 3 sources that are underestimated corresponds to small eigenvalues. Figure 6.10 compares the underestimation by 3 sources (denoted by CP-3) and 6 sources (denoted by CP-6). It is found that all the results are comparable, which means that CP is robust to underestimation of the number of sources.

Conclusion of the chapter

CP and FISTA have been compared from different perspectives. First, CP and FISTA can be understood in an uniform way from a signal reconstruction perspective. From the practical choice of an algorithm, it is suggested to use CP in cases with low SNRs and small numbers of sources, while FISTA is suggested in cases with high SNRs and large numbers of sources. In this chapter, it is also concluded that CP is robust to the underestimation and overestimation of the numbers of sources.

Conclusion and future works 7.1 Summary of the contribution

The problem of sequential measurements without reference has been investigated in this thesis. Compared with simultaneous measurement, the phase relationships between consecutive positions is missing and results in missing entries of the spectral matrix. Thus, the problem boils down to a spectral matrix completion problem of which a unique solution can exist only by providing additional constraints. Here, the realistic assumption of a sparse eigenvalue spectrum has been considered. Firstly, a structured low rank model has been formulated to find a full spectral matrix subject to reduced rank, measurements fitting, constraints of hermitian symmetry, and the spatial continuity of the sound field : which can be iteratively solved by CP as

S k = P C 1 P C 2 P C 3 (S k-1 ). (7.
2)

The CP algorithm has been proposed based on the idea of iterations between three predefined sets. The algorithm is easy to implement, yet it has no guarantee to find the global minimum and the number of acoustic sources needs to be known a priori or estimated. The rate of convergence of CP is linear. It has been found to perform very well in scenario with low SNRs and a relatively small number of sources.

Besides, a weakly sparse eigenvalue spectrum model was formulated which led to the minimization of the nuclear norm of the spectral matrix subject to data missing of spectral matrix, given constraints of hermitian symmetry and spatial continuity of the sound field : The FISTA has been proposed to solve the nuclear norm minimization, which is recognized as convex optimization with global minimum in guarantee. It has many parameters to tune and the number of sources is not necessary known. The rate of convergence of FISTA is O(1 k 2 ) which places it in the family of fast methods. It has been found that FISTA performs well in scenario with high SNRs and relatively large numbers of sources. CP and FISTA have been validated in different simulation scenarios : sequential measurements without and with references, with various SNRs, with different shift distances, and shown to achieve comparable performance as the MS method. Thus, CP and FISTA may be considered as alternative methods to MS in the configuration of with references. Both CP and FISTA are based on the following two physical assumptions :

1) Low rank spectral matrix. The high correlation of the acoustical field implies a low rank of the spectral matrix. The rank of the spectral matrix is consistent with the number of "virtual sources", i.e. the number of equivalent uncorrelated sources the acoustical field is composed of.

2) Continuity of the acoustic field. The high correlation of the acoustical field implies that two columns of the spectral matrix tend to become identical as the corresponding microphone positions become arbitrarily close to each other.

A big picture about the conclusion of thesis is depicted in Fig. 7.1. Finally, the three main results obtained in this thesis are highlighted hereafter :

• CP and FISTA can be both used to realize sequential measurements without and with references ;

• CP and FISTA can be considered as a better alternative to MS when the number of references is insufficient or when references are of poor quality ;

• CP is suggested to be used in cases with low SNRs and relatively small numbers of sources, while FISTA is suggested to be used in cases with high SNRs and relatively large numbers of sources.

Future works

Three main orientations will be considered in future works : Noise modeling, More complicated measurements setups, Model parameters optimization. Noise modeling : In the current work, measurement noise has not been explicitly modelled so that it results in estimation errors which are spread in the final recovered matrix S. In future work, measurement noise could be modelled as a diagonal matrix D and the weakly sparse eigenvalue spectrum model reformulated as follows. Generalize current methods to more complicated measurements setup : in this thesis, the sequential measurements have been carried out in the same plane. Assuming that the acoustic sources are located inside a 3D closed object, sequential measurements would rather be performed in a way to enclose the object of interest. This is illustrated in Fig. 7.2, where sequential measurements are implemented so as to cover completely a 3D object. In order to adapt the proposed approaches to such configurations, other spatial basis need to be investigated that are better suited to the 3D case. Model parameters optimization : In acoustical reconstruction, the sequential measurements was formulated as the physical model with different physical parameters, e.g. the topology of the source surface, the a priori spatial distribution of the sources, the geometry of the array, the relative shift distances between the sequential measurements etc. When the sequential measurements was formulated as the spectral matrix completion and solved by the proposed methods in this thesis, which is the data model. The relations between the optimal algorithm parameters and a given physical parameters setup are complex, and How to choose the optimal parameters in the algorithms for a given physical parameters setup will be an interesting topic in the future, which the idea is illustrated in Fig. Frobenius norm of matrix X ⟨x, y⟩ inner product diag(X) diagonal values of matrix X Rank(X) rank of matrix X trace(X) trace of the matrix X prox proximity operator shrink shrinkage operator ∇f (x) gradient of function f (x) X ⪰ 0 matrix X is positive semidefinite P C (x) = y y is attained by projecting x onto the set C.

1 Figure 1 . 1 -

 111 Figure 1.1 -(a) Measurement with a prototype array, (b) Sequential measurements for making the array denser, (c) Sequential measurements for making the array larger.

Figure 1 . 2 -

 12 Figure 1.2 -Spectral matrix Ŝm pp (left, data missing matrix), unknown full spectral matrix to be searched (right, full matrix).

  Figure 1.3 -A big picture about this thesis.

Figure 2 . 1 -

 21 Figure 2.1 -Illustration of the position of the approach proposed in this thesis with respect to the state-of-the-art.

Figure 2 . 2 -

 22 Figure 2.2 -Illustration of the conceptions of different acoustical imaging techniques : NAH, Patch NAH, Compressed sensing NAH.

Figure 2 . 3 -

 23 Figure 2.3 -Illustration of the assumptions underlying different acoustical imaging techniques : NAH, Patch NAH, Compressed sensing NAH.

Figure 2 . 4 -

 24 Figure 2.4 -Illustration the difference between intrusive and non-intrusive methods.

Figure 3 . 1 -

 31 Figure 3.1 -Sequential measurements versus simultaneous measurement (measurement mode).

Figure 3 .

 3 2 illustrate the conception how to generate the snapshots of sound pressure signal, see Ref[START_REF] Pintelon | System identification : a frequency domain approach[END_REF] [64][START_REF] Antoni | The spectral kurtosis : a useful tool for characterising non-stationary signals[END_REF]

Figure 3 . 2 -

 32 Figure 3.2 -illusion of the conception of snapshots that are calculated by Short Time Fourier Transform (STFT) tapered with a smooth data window.

Figure 3 . 3 -

 33 Figure 3.3 -left : Sequential measurements spectral matrix Ŝm pp without reference (data missing matrix), middle : Sequential measurements spectral matrix Ŝm pp with references (data missing matrix), right : Simultaneous measurement spectral matrix Ŝpp (full matrix).

Figure 3

 3 Figure 3.4 -(a) Full spectral matrix S L pp , (b) Eigenvalues distribution of full spectral matrix S L pp , (c) Measured spectral matrix Ŝm pp , (d) Eigenvalues distribution of measured spectral matrix Ŝm pp .

Algorithm 1 : 1 : 2 : 4 :σ 2 i u i u * i . 5 :

 11245 Cyclic Projection (CP)/(the Matlab codes are attached in the appendix section) Starts with Ŝ{0} ∈ C M P ×M P . While Iteration number < M max do 3: S{k} = P C 1 ( Ŝ{k} ). see Eq.(3.21) for calculation (Keep the measurements unchanged in diagonal blocks : project from set C 3 onto set C 1 ) S{k} = ∑ n i=1 Š{k} = P C 2 ( S{k} ) = ∑ r i=1 σ 2 i u i u * i . (Eigenvalues truncation : project from set C 1 onto set C 2 ) 6: If Stopping criteria ≤ SC, break 7: Ŝ{k+1} = P C 3 ( Š{k} ) = Ψ Š{k} Ψ * .

Figure 4 . 1 -

 41 Figure 4.1 -(a) Eigenvalues distribution of spectral matrix with SNR = 30 dB and 3 sources, (b) Eigenvalues distribution of spectral matrix with SNR = 30 dB and 25 sources. (the simulation setup in this example is explained in Fig. 5.1 (b).)

Figure 4 .

 4 1 (a) shows the eigenvalues distribution of the spectral matrix in the case of SNR = 30 dB and 3 sources, the rank r = 3 is clearly identifiable.

Figure 4 .

 4 1 (b) shows the eigenvalues distribution of the spectral matrix in the case of SNR = 30 dB and 25 sources, the determination of the rank is less obvious. In order to overcome this difficulty, the structured low rank model can be reformulated as a weakly sparse eigenvalue spectrum model as follows : minimize S S * subject to A (S) -Ŝm pp F ≤ ΨSΨ * -S F ≤ ε S * = S ⪰ 0

. 25 )

 25 Each projection P C (x) = y means that y is attained by projecting x ∈ R n onto the nonempty closed convex set C. It can be interpreted as the solution of the following optimization problem minimize y I C (y) + xy 2 2 , (4.26)

.46) 4 . 4 FISTAFigure 4 . 2 -A (S) -Ŝm pp 2 F

 44422 Figure 4.2 -Illustration the iterations of FISTA.

(4. 47 ) 2 F

 472 It is noted that ΨSΨ * -S F is a small quantity that reflects the representation error between ΨSΨ * and S ; thus this term can be simply enforced by setting ΨSΨ * = S in each iteration. This can be further written as minimize S A (S) -Ŝm pp + λ S * subject to x ∈ C3, (4.48) where C 3 = {B ∈ C M P ×M P ∶ B = ΨAΨ * , ∀A = A * }. The corresponding projection P C 3 (A) denotes an operator that maps the elements (A) onto set C 3 by B = ΨAΨ * . As compared with the generic model (4.36), it is found that A (S) -Ŝm pp 2

pp 2 F

 2 is given by ∇f (S k-1 ) = A (S) -Ŝm pp . (4.50)

Figure 5 . 1 -

 51 Figure 5.1 -Simulation setups : (a) A prototype array with 25 microphones (center microphone position is marked with red cross), (b) Setup 1, (c) Setup 2, (d) Setup 3 (the references microphones are marked with red point), (e) Setup 4, (f) Setup 5 (only the case without reference is shown in setups 3, 4, and 5, and center microphone position in prototype array is shifted 9 times sequentially which is marked with red crosses in setups 1, 2, 3, 4 and 5).

Figure 5 . 2 -

 52 Figure 5.2 -Velocity distribution of simulated sources.

  [START_REF] Johnson | Array signal processing : concepts and techniques[END_REF].1 (c)) at positions (-0.20, 0.25), (-0.20, 0.15), (-0.20, 0.05), (-0.20, -0.05).-3 references at positions (-0.20, 0.25), (-0.20, 0.15), (-0.20, 0.05).-2 references at positions (-0.20, 0.25), (-0.20, 0.15).-1 reference at positions (-0.20, 0.25) (meters).

  [START_REF] Johnson | Array signal processing : concepts and techniques[END_REF].1 (f)) : the prototype array is moved 9 times sequentially at positions (-0.22, -0.22), (-0.22, 0.08), (-0.22, 0.22), (0.08, -0.22), (0.08, 0.08), (0.08, 0.22), (0.22, -0.22), (0.22, 0.08), (0.22, 0.22) (meters).The parameters of the CP algorithm are chosen as follows : maximum iteration number M max = 6000, stopping criteria SC = 10 -3 , rank r = 3. The parameters of the FISTA are chosen as follows : maximum iteration for each regularization step N m = 10, stopping criteria SC = 10 -3 , initial regularization parameter λ 0 = Ŝm pp F , final regularization parameter λ d = 1e -16 λ 0 , ratio η = 0.7 and step size µ = 1.2. Two figures of merit are considered : (1) the spectral matrix completion error (MCE) calculated as MCE

Figure 5 . 3 -Figure 5 . 4 -

 5354 Figure 5.3 -Sequential measurements with different spatial bases : (a) Lx = Ly = 1 m, (b) Lx = Ly = 2 m (spatial bases are indicated by their dimensions 225 × 64, 225 × 81, 225 × 100, 225 × 121, 225 × 144), (c) 64 × 64 Fourier basis with parameters Lx = Ly = 1 m and ∆x = ∆y = 0.12 m, (d) 64 × 64 Fourier basis with parameters Lx = Ly = 2 m and ∆x = ∆y = 0.24 m.

Figure 5 . 5 -

 55 Figure 5.5 -Sequential measurements without reference (Setup 1) : (a) MCE without reference by CP, (b) ARE without reference by CP.

Figure 5 . 6 -

 56 Figure 5.6 -Sequential measurements without reference (Setup 1) : (a) MCE without reference by FISTA, (b) ARE without reference by FISTA.

Figure 5 .Figure 5 . 7 -

 557 Figure 5.7 -Reconstruction results analysis with CP : (a) MCE without reference, (b) ARE without reference. MCE and ARE : (c)(d) with 1 reference, (e)(f) with 2 references.

Figure 5 . 8 -Figure 5 . 9 -

 5859 Figure 5.8 -Reconstruction results analysis with CP, MCE and ARE : (a)(b) with 3 reference, (c)(d) with 4 references.

Figure 5 . 10 -

 510 Figure 5.10 -Reconstruction results analysis with FISTA, MCE and ARE : (a)(b) with 3 reference, (c)(d) with 4 references.

Figure 5 . 11 -Figure 5 . 12 -Figure 5 . 13 -Figure 5 . 14 -Figure 5 . 15 -

 511512513514515 Figure 5.11 -Reconstruction error with different SNRs and different numbers of references : (a) MCE with SNR = 60 dB, (b) ARE with SNR = 60 dB, (c) MCE with SNR = 30 dB, (d) ARE with SNR = 30 dB. (CP results are shown in black (as the legend in (b)) : with 1 reference case by upward-pointing triangle ; with 2 references case by downward-pointing triangle ; with 3 references case by diamond ; with 4 references case by circle. MS results are shown in gray (as the legend in (b)) : with 1 reference case by square ; with 2 references case by right-pointing triangle ; with 3 references case by cross ; with 4 references case by asterisk.)In this section, the simulation based on Setup 1 to 5 with references are firstly investigated. CP and FISTA are also investigated in the case of references under setup 2 and

  5.11 (a) and Fig. 5.13 (a)) ; when SNR = 30 dB, CP and FISTA are comparable with MS (see Fig. 5.11 (c) and Fig. 5.13 (c)) ; CP/FISTA outperforms MS when SNR ≤ 20 dB (see Fig. 5.12 (a), Fig. 5.12 (c), Fig. 5.14 (a) and Fig. 5.14 (c))

  5.11 (a) and Fig.5.13 (a)) ; CP/FISTA outperforms MS when the SNR ≤ 30 dB(see Fig.5.11 (c), Fig.5.12 (a), Fig.5.12 (c), Fig.5.13 (c), Fig.5.14 (a) and Fig.5.14 (c)). Similar conclusion holds as in the previous case, yet with increased sensitivity to noise.

Fig. 5 .

 5 11 (a), Fig. 5.11 (c), Fig. 5.12 (a), Fig. 5.12 (c), Fig. 5.13 (a), Fig. 5.13 (c), Fig. 5.14 (a) and Fig. 5.14 (c).

  5.11 (b), Fig. 5.11 (d), Fig. 5.12 (b), Fig. 5.12 (d), Fig. 5.13 (b), Fig. 5.13 (d), Fig. 5.14 (b) and Fig. 5.14 (d) is not only influenced by the spectral matrix completion errors but also the regularization. The spectral matrix completion errors are not always positively correlated with the acoustical source reconstruction error : ARE could be low with a high MCE due to the results of regularization.

Figure 5 .

 5 Figure 5.16 -(left) Experimental setup, (right) Spectral matrix completion error and acoustical source reconstruction error between CP/FISTA and MS methods.

  ,

  Figure 5.21 shows the spectral matrix completion results at different frequencies with different methods : (a) Spectral matrix completion result at 3000 Hz by CP, (b) FISTA and (c) MS ; (d) Spectral matrix completion result at 4000 Hz by CP, (e) FISTA and (f) MS ; (g) Spectral matrix completion result at 5000 Hz by CP, (h) FISTA and (i) MS ; (j) Spectral matrix completion result at 6000 Hz by CP, (k) FISTA and (l) MS. All these experimental results demonstratethat both CP and FISTA can be successfully applied to realize sequential measurements without references, and they can obtain very competitive results as compared with MS method in the setup of with references, which incurs the extra cost for the acquisition system. Thus, CP and FISTA may be considered as a better alternative to MS in practice.

Figure 5 .

 5 Figure 5.17 -(a)-(i) Acoustical source reconstruction at 3000 Hz for a fixed position of the array, (j) with the MS method, (k) the CP method, and (l) FISTA.

Figure 5 .

 5 Figure 5.18 -(a)-(i) Acoustical source reconstruction at 4000 Hz for a fixed position of the array, (j) with the MS method, (k) the CP method, and (l) FISTA.

Figure 5 .

 5 Figure 5.19 -(a)-(i) Acoustical source reconstruction at 5000 Hz for a fixed position of the array, (j) with the MS method, (k) the CP method, and (l) FISTA.

Figure 5 .

 5 Figure 5.20 -(a)-(i) Acoustical source reconstruction at 6000 Hz for a fixed position of the array, (j) with the MS method, (k) the CP method, and (l) FISTA.

Figure 5 .

 5 Figure 5.21 -(a) Spectral matrix completion result at 3000 Hz by CP, (b) FISTA and (c) MS ; (d) Spectral matrix completion result at 4000 Hz by CP, (e) FISTA and (f) MS ; (g) Spectral matrix completion result at 5000 Hz by CP, (h) FISTA and (i) MS ; (j) Spectral matrix completion result at 6000 Hz by CP, (k) FISTA and (l) MS.

Figure 6 . 1 -

 61 Figure 6.1 -Understanding CP/FISTA from a signal reconstruction perspective.

Figure 6 . 2 -

 62 Figure 6.2 -Comparison between CP and FISTA (eigenvalue spectrum perspective).

Figure 6 . 4 -Figure 6 . 5 -

 6465 Figure 6.4 -A prototype array with 25 microphones (shown in previous chapter) are shifted (a) 16 times, (b) 25 times, (c) 36 times sequentially which is marked with red crosses.

  Number of sources is 25 : when SNR = 60 dB and 30 dB, FISTA outperforms CP in the frequency range below 2000 Hz and CP outperforms FISTA in the frequency range above 2000 Hz. When SNR = 20 dB, FISTA is comparable with CP. When SNR = 10 dB, CP outperforms FISTA. • Number of sources is 36 : FISTA outperforms CP when SNR = 60 dB and SNR = 30 dB as seen in Fig. 6.6 (d) ; when SNR = 20 dB, FISTA is comparable with CP. When when SNR = 10 dB, CP outperforms FISTA.

Figure 6 . 6 -Figure 6 . 7 -Figure 6 . 8 -Figure 6 . 9 -Figure 6 . 10 -

 66676869610 Figure 6.6 -MCE with different SNRs and numbers of sources : (a) MCE with 3 sources, (b) MCE with 9 sources, (c) MCE with 16 sources, (d) MCE with 25 sources, (e) MCE with 25 sources. (CP results are shown in black (as the legend in (b)) : SNR = 60 dB ○ ; SNR = 30 dB ◇ ; SNR = 20 dB ▽ ; SNR = 10 dB △. FISTA results are shown in gray (as the legend in (b)) : SNR = 60 dB ○ ; SNR = 30 dB ◇ ; SNR = 20 dB ▽ ; SNR = 10 dB △.)

  find S ∈ C M P ×M P such that ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ Rank(S) = r (reduced rank) A (S) -Ŝm pp F ≤ (measurements fitting)ΨSΨ * -S F ≤ ε (spatial continuity of the sound field) S * = S ⪰ 0, (hermitian symmetry and positive semi-definiteness) (7.1)

  to A (S) -Ŝm pp F ≤ (measurements fitting) ΨSΨ * -S F ≤ ε (spatial continuity of the sound field) S * = S ⪰ 0, (hermitian symmetry and positive semi-definiteness) (7.3)which can be iteratively solved by FISTA as (acceleration steps are omitted in this conclusion)S k = P C 3 (shrink(S k-1 -µ∇f (S k-1 ), λµ)).(7.4)

Figure 7 . 1 -

 71 Figure 7.1 -Conclusion of the thesis.

S

  * + D l 1 subject to A (S) -Ŝm pp F ≤ (measurements fitting) ΨSΨ * -S F ≤ ε (spatial continuity of the sound field) S * = S ⪰ 0 (hermitian symmetry and positive semi-definiteness) d i,j = 0 if i ≠ j, (noise modeling) (7.5) with d i,j , i, j = 1, ...M P the elements of D. The noise modeling can be also considered in the Bayesian framework : by means of a Wishart likelihood function with structured covariance matrix and solved by Variational Bayes (VB) or Markov Chain Monte Carlo (MCMC) methods.

Figure 7 . 2 -

 72 Figure 7.2 -Generalize current methods to more complicated measurements setup.

7

 7

  .

Figure 7 . 3 -

 73 Figure 7.3 -Model parameters optimization.

  mu = max(mu * 0.7, mu_final); for i = 1:N_iter Gk = Yk + Li.* (D_datamissing -Yk.*PRDmatrix); % gradient descent Gk = psi_B*Gk*psiB_T; % smooth the spectral matrix % low rank estimation by SVD % Gk = (Gk + Gk')/2; % [U,S,V] = svd(Gk); % Xkk = U*diag(max(diag(S) -mu*Li,0))*V'; % low rank estimation by EVD Gk = (Gk + Gk')/2; [V,S]= eig(Gk); [Y,I]= sort(diag(S),'descend'); S = diag(Y); V = V(:,I); Xkk = V*diag(max(diag(S) -mu*Li,0))*V'; tk1 = 0.5 + 0.5*sqrt(1+4*tk^2); Yk = Xkk + ((tk -1)/tk1)*(Xkk -Xk); % low rank matrix update Xk = Xkk; tk = tk1; end err = norm(D_datamissing -Yk.*PRDmatrix,'fro')/norm(D_datamissing,'fro'); if err < SC; break end total_iter = N_iter * outer_iter; outer_iter = outer_iter + 1; end R_matrix = Xk; % completed spectral matrix para_val.cost_time = toc; % cost time para_val.mu = mu; % final mu para_val.err = err; % residual error para_val.total_iter = total_iter; % total iteration steps

Table 2 .

 2 2 -Conclusion of the current methods.

		proaches	Ap-	Assumption	Advantages / Dis-advantages
					Easy to implement /
	Single measure-ment	Patch NAH [27] [28] [29] [35] [33]	Band-limitation	Errors increase with the extrapolation sur-
					face
					Rooted on strong
	Single measure-ment	Compressed Sensing NAH [45] [44] [46] [47]	Sparsity	mathematical theory / Limited application due to the dictionary
					design
	Sequential mea-surements with references	with references meth-ods [56] [57] [58]	Stationary of the acoustical field	Very efficient / Incur extra cost
	Sequential mea-surements with-out reference	EM [62]		Stationary and field high correlation of the acoustical	Very efficient / Intru-sive method

End for if 13

  : λ k = max(ηλ k-1 , λ d ) 14: Go to step 2.

	1 )
	9: End
	10: If Stopping criteria ≤ SC,
	11: break
	12:

Table 6 .

 6 1 -Comparison between CP and FISTA.

	CP	FISTA

Cette thèse est accessible à l'adresse : http://theses.insa-lyon.fr/publication/2015ISAL0023/these.pdf © [L. Yu], [2015], INSA de Lyon, tous droits réservés

Acknowledgments

Eigenvalues Filter 2 r truncation

Simulation and parameter analysis

Convergence and computational cost

The numerical convergence and computation cost of CP and FISTA are investigated in this section. The simulation is carried out with Setup 1 and the measurements are simulated without noise at a given frequency 1200 Hz, using the workstation with an Intel core i3 cpu 530 (2.93Ghz,4G memory), Windows 7 and Matlab R2013b. The parameters of the CP algorithm are chosen as follows : maximum iteration number M max = 1680, stopping criteria SC = 10 -4 , rank r = 3. The parameters of the FISTA are chosen as follows : maximum iteration for each regularization step N m = 10, stopping criteria SC = 10 -4 , initial regularization parameter λ 0 = Ŝm pp F , final regularization parameter λ d = 1e -26 λ 0 , ratio η = 0.7 and step size µ = 1.2 (the maximum iteration number M max = 1680 for the given parameters). Figure 6.5 (a) shows a numerical example that the MCE decreases with the number of iteration, and the MCE of FISTA iteration decreases faster than the CP. Figure 6.5 (b) shows the computational time ([s] : second unit) increases with the number of sequential measurements from 9, to 16, to 25, to 36, and the corresponding spectral matrix size to complete is 225 × 225, 400 × 400, 625 × 625, 900 × 900. Figure 6.5 (c) shows the MCE for CP and FISTA with the increase of the number of sequential measurements. In Figs. [START_REF] Leclere | Acoustic imaging using under-determined inverse approaches : Frequency limitations and optimal regularization[END_REF].5(b)(c), the CP and FISTA use up all the 1680 iterations (with the same number of iteration used), it is found that FISTA can achieve better MCE and expend less computational times with the number of sequential measurements increases. Une limitation fondamentale du problème inverse acoustique est déterminée par la taille et la densité de l'antenne de microphones. Une solution pour atteindre une grande antenne et / ou à forte densité de microphones est de scanner l'objet d'intérêt en déplaçant séquentiellement une antenne de petites dimensions (mesures dites séquentielles). La différence entre des mesures séquentielle et une mesure simultanée en tous points est que dans ce dernier cas, l'intégralité de la matrice interpectrale peut être estimée, contrairement au cas des mesures séquentielles qui ne permettent l'estimation que d'une partie réduite de cette matrice; les éléments croisés (interspectres) entre deux points de mesures n'appartenant pas à la même séquence ne sont pas estimés. Néanmoins, ces données restent nécessaires pour la reconstruction acoustique. Dans l'approche classique, une ou plusieurs références sont utilisées pour retrouver les données manquantes. L'objet de cette thèse est de récupérer les éléments manquants de la matrice interspectrale sans capteurs de référence, dans le cas où le champ acoustique est suffisamment cohérent pour mettre en oeuvre les mesures séquentielles. Deux modèles de spectre de valeurs propres parcimonieux sont proposés pour résoudre ce problème, le premier impose une valeur faible de rang, tandis que le second repose sur la minimisation de la norme nucléaire (recherche d'une solution faiblement parcimonieuse).

Matlab codes of CP and FISTA

MOTS-CLÉS :

reconstruction inverse acoustique, mesures séquentielles sans référence, modèles de spectre de valeurs propres parcimonieux, projection cyclique