# Modélisation, caractérisation et exploitation des milieux réverbérants

Florian Monsef

► To cite this version:

Florian Monsef. Modélisation, caractérisation et exploitation des milieux réverbérants. Electromagnetism. Université Paris XI Paris-Sud, 2019. tel-04552730

HAL Id: tel-04552730

https://hal.science/tel-04552730

Submitted on 19 Apr 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

Public Domain
NNT : 2016SACLCxxx

HABILITATION À DIRIGER DES RECHERCHES
de

l'Université PARIS-SACLAY
Laboratoire Geeps

(Univ. Paris-Sud, CNRS, CentraleSupélec, UPMC Sorbonne Université)

Modélisation, caractérisation et exploitation des milieux réverbérants

ÉCOLE DOCTORALE N° : 575

Physique et Ingénierie : Électrons, Photons, Sciences du vivant (EOBE)

Spécialité : Electromagnétisme

Par

Florian MONSEF

Habilitation à diriger les Recherches

Habilitation présentée et soutenue à Gif-sur-Yvette, le 17 juin 2019

## Composition du Jury :

### Rapporteurs :

Julien De Rosny Directeur de Recherche CNRS, Institut Langevin, ESPCI

Philippe Besnier Directeur de Recherche CNRS, IETR INSA Rennes

Elodie Richalot Professeur des Universités, Univ. Paris-Est Marne-la-vallée

### Examinateurs :

Marc Lambert Chargé de Recherche CNRS, Geeps, Univ. Paris-Sud, CentraleSupelec

Fouzia Boukour Directrice de Recherche IFSTTAR, Villeneuve d'Ascq

Franco Moglie Professeur à l'Université des Marche, Ancona (Italie)

Cyril Decroze Maître de Conférences à l'Université de Limoges (XLIM)
<figure><img src="image_1.png" /></figure>
Table des matières

<table><thead><tr><th></th><th></th><th></th></tr></thead><tbody><tr><td>1</td><td>Curriculum Vitae</td><td>1</td></tr><tr><td>I</td><td>Formation universitaire</td><td>1</td></tr><tr><td>II</td><td>Parcours Professionnel</td><td>1</td></tr><tr><td>III</td><td>Domaines de Compétences</td><td>1</td></tr><tr><td>III.1</td><td>Compatibilité Electromagnétique (CEM)</td><td>2</td></tr><tr><td>III.2</td><td>Acoustique des Chambres</td><td>2</td></tr><tr><td>III.3</td><td>Retournement Temporel (RT)</td><td>2</td></tr><tr><td>IV</td><td>Activités d'enseignement</td><td>2</td></tr><tr><td>V</td><td>Activités d'encadrement</td><td>3</td></tr><tr><td>V.1</td><td>Stagiaires encadrés</td><td>4</td></tr><tr><td>Stagiaires encadrés/co-encadrés</td><td></td><td>4</td></tr><tr><td>V.2</td><td>Doctorants encadrés</td><td>4</td></tr><tr><td>VI</td><td>Collaborations</td><td>4</td></tr><tr><td>2</td><td>Physique des Chambres Réverbérantes</td><td>5</td></tr><tr><td>I</td><td>Régime surmodé ou régime diffus</td><td>6</td></tr><tr><td>II</td><td>Régime non surmodé</td><td>8</td></tr><tr><td>II.1</td><td>Ergodicité</td><td>8</td></tr><tr><td>II.2</td><td>Approche modale</td><td>11</td></tr><tr><td>III</td><td>Du Régime sous-modé au régime sur-modé</td><td>19</td></tr><tr><td>IV</td><td>Conclusion et perspectives</td><td>22</td></tr><tr><td></td><td>Références</td><td>22</td></tr><tr><td>3</td><td>Retournement Temporel</td><td>65</td></tr><tr><td>I</td><td>Retournement Temporel classique</td><td>65</td></tr><tr><td>II</td><td>Métrique utile : le contraste en énergie et le contraste pic</td><td>68</td></tr><tr><td>III</td><td>Contributions au RT classique</td><td>69</td></tr><tr><td>III.1</td><td>Influence du nombre de sources</td><td>69</td></tr><tr><td>III.2</td><td>RT pour la synthèse de front d'ondes : la TREC</td><td>73</td></tr><tr><td>IV</td><td>Conclusion et perspectives</td><td>80</td></tr><tr><td></td><td>Références</td><td>80</td></tr><tr><td>4</td><td>Mesure d'antenne en milieu échogène</td><td>109</td></tr><tr><td>I</td><td>Modélisation du problème dans le cas monoécho</td><td>109</td></tr><tr><td>II</td><td>Pertinence du modèle proposé</td><td>111</td></tr><tr><td>II.1</td><td>Hypothèse de Champ lointain</td><td>111</td></tr><tr><td>II.2</td><td>Diffraction de l'antenne</td><td>113</td></tr></tbody></table>
Table des matières

<table><tbody><tr><td></td><td>II.3 Plaque de dimension finie</td><td>115</td></tr><tr><td>III</td><td>Inversion du problème</td><td>117</td></tr><tr><td></td><td>III.1 Nombre de conditionnement dans l'approche classique</td><td>117</td></tr><tr><td></td><td>III.2 Approche par Régresseurs</td><td>119</td></tr><tr><td>IV</td><td>Configurations multi-échos</td><td>124</td></tr><tr><td></td><td>IV.1 Structures Diédrales</td><td>125</td></tr><tr><td></td><td>IV.2 Cas des plaques parallèles</td><td>125</td></tr><tr><td>V</td><td>Conclusion et perspectives</td><td>126</td></tr><tr><td></td><td>Références</td><td>126</td></tr><tr><td>5</td><td>Liste de mes publications</td><td>127</td></tr><tr><td>I</td><td>Articles dans des revues internationales à comité de lecture</td><td>127</td></tr><tr><td>II</td><td>Articles dans des conférences internationales à comité de lecture avec actes</td><td>128</td></tr><tr><td>III</td><td>Brevets</td><td>129</td></tr></tbody></table>
2.1 Chambre réverbérante à brassage de modes de CentraleSupelec. Le pale en mouvement modifie les conditions limites et permet de « brasser le champ » au sein de l'enceinte. 6

2.2 Taux de réjection d'hypothèse de gaussianité d'une distribution à N échantillons corrélés. Le niveau de corrélation (type lag-1) est représenté en abscisse, le taux de réjection est donné par l'axe des ordonnées. Les tests d'hypothèse suivants ont été étudiés : Liliefors, Anderson-Darling (AD), Kolmogorov-Smirnov (KS) ainsi que le $\chi^2$. 8

2.3 Chambre réverbérante acoustique du LNE (Trappes). Les plaques suspendues sont des diffuseurs servant à casser tout degré de symétrie afin d'assurer une chaoticité du milieu; les microphones suspendus servent quant à eux à échantillonner le champ. Au centre une tige rotative qui sert à échantillonner le champ; le mouvement de cette tige a un effet de brassage du champ. 9

2.4 (a) Illustration du recouvrement modal (bas à gauche) où chaque mode est représenté. Réponse spectrale résultante (gauche haut). (b) Réponse expérimentale montrant la complexité effective de la superposition modale. $B_e$ est la bande effective à considérer incluant le nombre utile de modes. 12

2.5 Valeur moyenne (trait gras) et limites des intervalles de confiance de $\hat{\zeta}_P^2$ fonction du nombre moyen $M_M$ de modes recouverts pour différentes valeurs de $N_{sp}$ positions indépendantes de brasseur (reportées sur le graphe). Les limites supérieures et inférieures de l'intervalle de confiance sont au-dessus et en-dessous du trait gras respectivement. 13

2.6 Densité de probabilité de Wigner (trait gras plein), de Poisson (tirets), des plus proches voisins ($NN$) (trait plein fin) et triangulaire (traits en pointillés). 15

2.7 Variabilité estimée $\hat{\zeta}_W^2$ (axe y gauche) de la densité d'énergie électrique, fonction du nombre $M_M$ de modes dans le cas où la chambre est vide (a) et chargée (b). Les déviations relatives à la valeur asymptotique de $1/3$ sont reportées (axe des ordonnées de droite). Les résultat expérimentaux (symboles gris), ceux de simulation MC (points noirs) et les résultats analytiques sont reportés (trait plein). Les barres d'incertitude verticales traduisent un intervalle de confiance de 95% et les barres horizontales rendent compte de la dispersion que le paramètre $M_M$ autour de sa valeur moyenne. 16

2.8 (a) Temps de réverbération estimés lorsque la chambre est vide (cercles) et chargée (carrés) (gauche). (b) Estimation du nombre de modes recouverts dans une bande modale dans les cas respectifs où la chambre réverbérante est vide (cercles) et chargée (carrés) (droite). 17

2.9 17
3.8 Réponses impulsionnelles des différents milieux étudiés : en-dehors de la chambre (en haut à gauche), dans un couloir (en haut) droite), dans la chambre réverbérante (en bas à gauche) et dans le hall illustré Fig. 3.7. 72

3.9 Illustration du SNR obtenu par RT avec une fréquence de 5 GHz : rapport des deux pics considérés (écho et impulsion utile). 73

3.10 SNR obtenus (gauche) dans les différents milieux indoor étudiés et le gain qu'apporte le RT (droite) par rapport une transmission classique. 73

3.11 L'arche constituant le coeur de la TREC consiste en un échantillonneur spatial du champ suivant les coordonnées polaires θ et φ (a) - Un zoom met en évidence le système de courroie utilisé ainsi que le chariot accueillant la sonde électro-optique (b). 74

3.12 Topographies 3D (gauche) et coupe 2D (droite) obtenues avec la TREC. 75

3.13 Topographie de champ visé où la concentration de l'onde est confinée spatialement et localement plane. 76

3.14 Pour une largeur W=5 et un angle Θ = 30°, on illustre les fonctions propres résultantes avec les valeurs propres λₙ correspondantes. 77

3.15 Topographies des fonctions Slepian pour W = 5 et un angle Θ = 20°. 77

3.16 Topographie de champ obtenue en incluant 4 Slepians. 78

3.17 Evolution du rapport des parties cohérentes et fluctuantes en fonction du degré de perturbation illustré par le coefficient (1 - ρ²). Les résultats expérimentaux sont représentés par les symboles pour trois statistiques différentes de chambre; les modèles théoriques sont représentés par les courbes en tirets. 79

3.18 Evolution relative du constraste en fonction du niveau de perturbation de la chambre, et ce, pour trois statistiques différentes. 79

4.1 Antennes utilisés pour l'analyse du concept proposé : (a) dimensions des AST (b) Coefficients de réflection S₁₁ correspondants et bandes de fonctionnement. 110

4.2 Principe de la mesure exploitant les échos. (a) La sonde est placée à une distance r₀ de l'AST et l'écho parcourt une distance r₁. (b) Schéma équivalent utilisant le principe des images mettant en évidence la mesure de deux points du diagramme en une mesure. 111

4.3 Erreur ϵ<sub>FF</sub> en module (dB) pour le cas TM (isotrope) - (a) cas du dipole (b) cas de l'antenne biconique. 112

4.4 Evolution des erreurs sur le module (a) ϵ<sub>∥E<sub>z</sub>∥</sub> et la phase (b) ϵ<sub>∠E<sub>z</sub></sub>, pour le cas TM isotrope. 113

4.5 Erreur obtenue en espace libre [dB] sachant l'hypothèse de champ lointain - cas TE : (a) ϵ(E<sub>x</sub>)<sub>HWD</sub> (b) ϵ(E<sub>y</sub>)<sub>HWD</sub> (c) ϵ(E<sub>x</sub>)<sub>biconique</sub> (d) ϵ(E<sub>y</sub>)<sub>biconique</sub>. 114

4.6 Cartographie de l'erreur liée au champ diffracté par l'antenne sous test en mode TM; plaque en xₚ = 10λ₀ : (a) ϵ₀<sup>d</sup>, (dipôle) (b) ϵ₁<sup>d</sup> (dipôle) (c) ϵ₀<sup>d</sup> (biconique) (d) ϵ₁<sup>d</sup> (biconique). 116

4.7 Configuration considérée sachant une plaque de dimension finie afin d'en évaluer l'effet de diffraction par les bords. 117
# CHAPITRE 1

## Curriculum Vitae

### I Formation universitaire

1999-2002 **Thèse de doctorat**, Institut d'Electronique Fondamentale, Université Paris-Sud
Spécialité : Electronique

1998-1999 **D.E.A** en Electronique et Conception de Circuits Intégrés. Option micro-ondes.

1995-1998 Ecole Normale Supérieure de Cachan - Intégration en 3ème année.
Agrégation en Ingénierie Electrique (option : électronique).

### II Parcours Professionnel

1999-2002 Maître de Conférences au département de Mesures Physiques à l'IUT d'Orsay,
Université Paris-Sud

2002-2009 Professeur Agrégé à l'IUT de Ville d'Avray au Département de Génie Electrique
Rang : 1er.

1995-1998 Ecole Normale Supérieure de Cachan - Intégration en 3ème année.
Agrégation en Ingénierie Electrique (option : électronique).

### III Domaines de Compétences

Depuis mon doctorat, j'ai pu développé des compétences en matière d'enseignement, de développement et de recherche. Les compétences en enseignement ont pu se faire grâce à un poste en tant qu'enseignant à plein temps au niveau IUT mais également avec des publics de techniciens et ingénieurs venus de l'industrie en formation continue. Au niveau de ces enseignements, des compétences techniques ont été mises en oeuvre et développées dans le cadre de mise en place de travaux pratiques et de projets.

Au niveau de la recherche, mon poste MCF m'a amené à une conversion thématique passant de la physique des semiconducteurs à l'électromagnétisme dans le domaine des microondes dans des cavités chaotiques à taille humaine, appelées *Chambres Réverberantes à Brassage de Modes (CRBM)*. Ces dispositifs sont essentiellement utilisés dans les domaines de la compatibilité électromagnétique.

L'investissement dans ces thématiques m'a conduit également en début de poste vers l'acoustique où des problématiques analogues existent. En outre, la technique de retournement tem-
# 1 Curriculum Vitae

porel (RT) appliquée au cas des milieux complexes a permis de développer de nouvelles perspectives métrologiques pour les mesures en CRBM.

La dernière thématique abordée dans le cadre d'une thèse a porté sur la mesure de diagramme de rayonnement d'antenne. Cette thématique est en semi-rupture avec les autres thématiques. Les recherches dans le cadre de cette thèse ont débouché également sur un dépôt de brevet et m'a permis d'appréhender de plus près les difficultés expérimentales liées à une « bonne » mesure d'antenne.

La liste qui suit donne une vision d'ensemble des thématiques abordées dans le domaine de la recherche et constituent le noyau central de mes compétences de chercheur.

## III.1 Compatibilité Electromagnétique (CEM)

Une partie de mes travaux ont porté sur la pertinence de certains postulats ou directives présentes dans les standards liés aux CRBM. Cela a permis de se familiariser avec les techniques propres à la CEM dans ce type d'environnement d'essais standardisés. Une fois les attentes de la communauté connues, une de mes motivations principale a porté sur une meilleure compréhension de la physique des CRBM dans un but d'affinage métrologique. Les compétences actuelles portent sur la capacité à prévoir les intervalles de confiance pour les mesures de puissance rayonnée.

## III.2 Acoustique des Chambres

Une autre partie de mes travaux a porté sur les problématiques abordées dans le cadre de chambre réverbérante acoustique. Une collaboration avec le Laboratoire National d'Essais (LNE) a permis de se familiariser avec un environnement standard et les protocoles de mesure dans ce domaine. L'immersion dans ce domaine permet de s'ouvrir à d'autres approches autour de problématiques analogues et permet un partage inter communautés. L'intérêt pour les travaux d'acoustique a permis un éclairage sur l'utilité d'introduire un niveau de chaoticité dans les chambres réverbérantes.

## III.3 Retournement Temporel (RT)

Familiarisation avec le RT classique développé en acoustique et prise en main du concept de RT généralisé développé au laboratoire. Cette technique est une façon de faire de la synthèse de front d'onde dans un environnement naturellement peu propice à ce type de synthèse. Une étude amont sur de la viabilité de cette technique pour de la mesure d'antenne et pour la génération d'impulsions de puissance en chambre réverbérante a été menée à travers des thèses de doctorat.

Au-delà de ces applications, le RT généralisé permet d'élargir la précision des mesures pour la CEM. Cet aspect mérite de plus amples investigations en cours.

# IV Activités d'enseignement

Mes activités d'enseignement ont essentiellement eu lieu au Département Mesures Physiques de l'IUT d'Orsay.
Dès mon recrutement, j'ai eu l'opportunité de monter un module portant sur la conversion photovoltaïque. Ce module a été chronophage car aucune expérience personnelle et locale n'était présente. Cela a consisté à mettre en place des procédures expérimentales de caractérisation de panneaux solaires d'une part et d'étudier le couplage de tels dispositifs avec un convertisseur de puissance de type hacheur. Il s'agit d'un module de 30h effectué en semestre 4 classique et en semestre décalé.

J'ai eu également la charge du module d'électronique d'instrumentation. Ma contribution dans le cadre de ce module a consisté à introduire des séances de travaux pratiques dédiées à la mesure et la mise en évidence du bruit thermique. J'ai eu également la responsabilité de modules d'automatique et de micro-contrôleur au début de mon poste. Certains sujets ont été remaniés. Il s'agit de modules de 30h effectué en semestre 3 et 4 ainsi qu'en semestres décalés.

Mon expérience d'enseignement s'est élargie aux écoles d'ingénieurs telles que l'Ecole Supérieure d'Electricité (aujourd'hui CentraleSupelec) et l'Ecole Supérieure de Physique et Chimie Industriels (ESPCI) de Paris.

Les périodes importantes jalonnant cette expérience d'enseignement sont résumées ci-dessous :

<table><thead><tr><td>Période</td><td>Nombre d'heures</td><td>Statut</td></tr></thead><tbody><tr><td>2009-2017</td><td>200h</td><td>Maître de Conférences IUT d'Orsay (Dpt Mesures Physiques), Université Paris-Sud. Domaines : automatique, électronique d'instrumentation et conversion photovoltaïque.</td></tr><tr><td>2002-2009</td><td>500h</td><td>Professeur Agrégé à l'IUT de Ville d'Avray au Département de Génie Electrique (électronique analogique, numérique et programmation).</td></tr><tr><td>1999-2002</td><td>64h</td><td>Moniteur à l'Université Paris-Sud en électronique</td></tr></tbody></table>

## V Activités d'encadrement

Dès mon arrivée à CentraleSupelec, j'ai été amené à participer à l'encadrement d'Houmam Moussa en profitant de l'immersion dans son sujet sur le RT pour démarrer ma conversion thématique. Mon introduction à l'encadrement de ce doctorant a eu lieu en 2009. Mon intégration s'est faite en rejoignant les thématiques développées par Andrea Cozza. Ces thématiques portent sur la physique des chambres réverbérantes et le retournement temporel (RT).

La deuxième année j'ai été amené à encadrer des stagiaires niveau Master ou équivalent. L'encadrement de doctorants a débuté avec le travail de Philippe Méton qui s'est déroulé officiellement avec la collaboration du CEA de Gramat. Cet encadrement s'est poursuivi avec le travail d'Henri Vallon en collaboration avec Thalès. Enfin, la thèse de Mouad Djedidi a été la dernière thèse encadrée.

Ci-dessous la liste des stagiaires et doctorants encadrés :
## V.1 Stagiaires encadrés

2012 Rana SALEH et Zareh DEMIRDJIAN, Université Libanaise, stage de fin d'études de 4 mois. Travail sur le retournement temporel.
Yann GIBONI, DUT de l'IUT de Cachan : échantillonneur linéaire du champ électronique.

2013 Djamila ABERKANE - Etudiante du Master CAT (Univ. Paris-Sud) - sujet : déconvolution des réponses en fréquence d'une CRBM pour en extraire la statistique des fréquences propres et des poids modaux

2014 Aloÿs DUFOUR, Stagiaire de Licence/Magistère de Physique Fondamentale. Etude de l'effet d'absorption de métamatériaux en CRBM.

## V.2 Doctorants encadrés

<table><thead><tr><td>Année</td><td>Taux Encadrement</td><td>Sujet et doctorant</td><td>Chapitre</td></tr></thead><tbody><tr><td>2011-2015</td><td>50%</td><td>Philippe METON, sujet : Étude et développement d'un concept de caractérisation rapide d'antennes basé sur le principe du retournement temporel en chambre réverbérante.</td><td>3</td></tr><tr><td>2012-2016</td><td>40%</td><td>Henri VALLON - Sujet : Focalisation d'impulsions forte puissance en utilisant le Retournement Temporel. Thèse CIFRE (Thalès Defense).</td><td>3</td></tr><tr><td>2013-2016</td><td>80%</td><td>Mouad DJEDIDI, sujet : Mesure de diagramme de rayonnement en milieu échogène l'effet d'absorption de métamatériaux en CRBM.</td><td>4</td></tr></tbody></table>

## VI Collaborations

<table><thead><tr><td>Collaborateur</td><td>Pays</td><td>Thématique</td></tr></thead><tbody><tr><td>Laboratoire National d'Essais (Garches)</td><td>France</td><td>Chambre réverbérante acoustique</td></tr><tr><td>Ramiro Serra (Université Eindhoven)</td><td>Pays-Bas</td><td>Test d'hypothèse en chambre réverbérante</td></tr><tr><td>Franco Moglie (Université d'Ancona)</td><td>Italie</td><td>Simulation chambre réverbérante</td></tr><tr><td>Thalès Defense</td><td>France</td><td>Mesures en chambre réverbérante en aluminium et financement thèse</td></tr><tr><td>CEA Gramat</td><td>France</td><td>Financement thèse</td></tr></tbody></table>
# CHAPITRE 2

## Physique des Chambres Réverbérantes

Les chambres réverbérantes à brassage de modes (CRBM) sont des dispositifs utilisés dans l'industrie, automobile et aéronautique notamment, pour des mesures de puissance rayonnée et des tests de susceptibilité électromagnétique; ce type d'applications s'inscrit dans les thématiques de la communauté de la Compatibilité Electromagnétique (CEM). La notion de brassage de modes est liée à l'utilisation d'un brasseur qui correspond à une pale métallique (rotative en général). Le rôle du brasseur est de générer un ensemble de scénarios de tests en modifiant les conditions limites.

Les CRBM sont des cavités microondes dont le volume peut aller de moins d'un mètre cube à quelques centaines de mètres cube. Un cliché de ce type de dispositifs est donné Fig. 2.1. Une CRBM est classée dans la catégorie des systèmes dits complexes car la prédiction du champ électromagnétique (EM) en leur sein y est délicate dans la pratique. En effet la position du brasseur ainsi que l'objet (forme, matière,...) sous test vont modifier et fixer la topographie finale du champ EM. Dans la pratique une centaine de positions de brasseurs sont utilisées typiquement. Ainsi, si l'on désire déterminer les topographies de champ par une approche déterministe l'utilisation de logiciels de simulation numérique est indispensable. La simulation de tels systèmes peut être néanmoins coûteux en temps de calcul.

Dans le cadre des travaux développés depuis 2009 le recours à la simulation pour des CRBM a fait l'objet de collaboration avec des numéritiens de l'université d'Ancona (Italie) dans le cadre de la thèse de Philippe Méton (cf. Chap. 3). Dans le reste des cas, une approche statistique du champ au sein de la CRBM a été privilégiée. Les avantages d'une approche statistique permet de prédire le niveau moyen du champ au sein de la CRBM sans se préoccuper des détails des divers éléments présents à l'intérieur de la chambre. La contrepartie est d'accepter un degré d'incertitude (intervalle de confiance) du niveau de champ en un point par rapport au niveau moyen prédit.

L'approche statistique a été initialement appliquée à un régime asymptotique pour lequel le nombre de modes excités dans la cavité est tel que le champ peut être vu en première approximation comme Gaussien. Nous présenterons notre contribution quant aux questions pendantes dans le cadre de ce régime et les questions soulevées par rapport aux approches existantes. Certaines références à l'acoustique y seront présentées dans la mesure où deux contributions ont été faites dans ce domaine. Dans un second temps nous abordons le degré d'homogénéité de la densité d'énergie et/ou de puissance au sein de la cavité.
# I Régime surmodé ou régime diffus

En utilisant une approche modale, le champ peut être décomposé en une superposition de modes. Il est d'usage en CEM de travailler principalement en régime harmonique, i.e., en appliquant une excitation sinusoïdale. Suivant la fréquence de travail, un ensemble plus ou moins important de modes seront excités, et ce, en vertu de la densité modale établie par Weyl (WEYL 1912). Cette dernière varie avec le cube de la fréquence pour une cavité 3D et est deux fois moins important dans le cas d'un champ de pression acoustique.

Dans une approche statistique, le champ est une somme de vecteurs aléatoires. Chacun de ces vecteurs correspond à un mode donné en un point considéré. Les excitations de ces modes sont fonction de la position de l'émetteur et du récepteur ainsi que de leur polarisation. Si le nombre d'éléments « significatifs<fnref n="1" /> » composant cette somme est suffisamment élevé, les parties réelles et imaginaires du champ, sur une composante cartésienne, pourront être vues comme des variables gaussiennes en vertu du théorème de la limite centrale. Si chaque mode est décomposé comme une somme d'ondes planes, alors le champ est classiquement représenté par une somme infinie d'ondes planes (HILL 1998). Cette représentation en acoustique est bien antérieure puisqu'elle était déjà proposé par Morse et Bolt dans les années 1940 (MORSE et BOLT 1944). Nous verrons un peu plus loin notre proposition, en rupture avec le modèle classique, de considérer un nombre fini d'ondes planes pour décrire le champ.

<figure><img src="image_2.png" /><figcaption>FIGURE 2.1 – Chambre réverbérante à brassage de modes de CentraleSupelec. Le pale en mouvement modifie les conditions limites et permet de « brasser le champ » au sein de l'enceinte.</figcaption></figure>

En électromagnétisme, on parle de champ gaussien, bien que la quantité concernée soit la partie réelle/imaginaire d'une composante cartésienne du champ. En acoustique en revanche, la grandeur d'intérêt correspond directement au champ de pression qui est une grandeur scalaire;

<footnote n="1">1. Nous verrons dans la section II.2 de ce que l'on entend par « significatifs »</footnote>
dans ce cas le champ de pression peut être effectivement vu comme gaussien. Les acousticiens parlent depuis les années 1950 de régime *diffus* (COOK, WATERHOUSE et al. 1955); la communauté CEM parle de régime *surmodé*.

Ainsi le modèle du champ gaussien reste une approximation asymptotique. Dans la pratique, une trentaine de termes suffise à adopter une approximation gaussienne ; ce nombre peut être vu à la baisse suivant la densité de probabilité des poids d'excitation. Lorsque ce régime est atteint à une fréquence donnée, la gaussianité du champ est, en première approximation<fnref n="2" />, assurée pour toutes les fréquences supérieures à celle-ci. Dans le domaine de l'acoustique cette fréquence est appelée fréquence de Schroëder - elle correspond à une excitation de trois modes dans la bande modale moyenne à -3 dB (SCHROEDER et K. KUTTRUFF 1962).

La contribution dans la recherche liée aux CRBM en régime diffus a porté sur l'utilisation d'outils statistiques permettant de vérifier la gaussianité du champ. Plus précisément, nous avons étudié la nécessité d'avoir des échantillons indépendants dans l'utilisation de tests d'adéquation.

Pour situer l'intérêt de cette étude, il convient de rappeler les deux points suivants. Premièrement, dans l'utilisation des CRBM, il est d'usage de rassembler les valeurs du champ (ou d'une de ses composantes) obtenues sur l'ensemble des positions du brasseur. Deuxièmement, il est de coutume de vérifier les performances du brasseur en estimant le niveau de corrélation entre chaque position du brasseur ; plus le niveau de corrélation est bas, plus le brasseur est vu comme performant dans sa capacité à générer des nouvelles topographies de champ. En-deça d'un certain niveau de corrélation, il est accepté<fnref n="3" /> de voir les échantillons comme décorrélés (et donc indépendants dans le cas gaussien). .

Il est souvent stipulé que les tests d'adéquation doivent être utilisés sur la base d'échantillons indépendants. Si cette règle s'avère vraie, cela poserait deux inconvénients quant à leur utilisation. D'une part, nous ne pourrions pas toujours être sûrs que lesdits échantillons sont vraiment indépendants (sachant l'incertitude liée à l'estimateur de corrélation) ; et d'autre part, cela pousserait à décimer les échantillons à disposition, augmentant la barre d'incertitude quant à l'issue du test d'adéquation. Dans le cas d'un test type $χ^2$ par exemple, il ne nous était pas intuitif que les échantillons doivent être indépendants sachant que le prérequis est que seules les observations, elles, soient indépendantes. Après quelques recherches, il nous est apparu qu'aucune preuve n'existe quant à la nécessité de satisfaire cette condition dans l'hypothèse nulle d'un test d'adéquation.

En collaboration avec Ramiro Serra de l'Université d'Eindoven (Pays-Bas), nous avons donc étudié (MONSEF, SERRA et al. 2015) l'influence de la corrélation sur l'issue des tests d'adéquation suivants : le $χ^2$, le test de Kolmogorov-Smirnov (KS), le test de le Liliefors et le test de d'Anderson-Darling (AD). L'étude a fait ressortir que les tests d'adéquation type Liliefors, AD et $χ^2$ sont insensibles à la corrélation entre échantillons jusqu'à des valeurs de 90% ! Le test KS est quant à lui très sensible à la corrélation. Pour les détails des travaux nous renvoyons à la Ref. (MONSEF, SERRA et al. 2015) donné en annexe. Toutefois, il est remarquable de constater la différence de comportement entre les tests KS et Liliefors. En effet, le test de Liliefors est un test KS qui prend en compte la moyenne et la variance de l'échantillon car ce test a été introduit dans le cas où le test d'adéquation s'adresse à un ensemble de données dont la moyenne et l'écart-type ne sont pas connus a priori. Le résultat principal est donné Fig. 2.2.

<footnote n="2">2. la densité modale n'est pas strictement monotone et des zones de déplétion modale peuvent exister comme clairement expliqué dans (COZZA 2012).</footnote>

<footnote n="3">3. d'après la norme IEC 61000-4-21 (Reverberation Chamber Test Methods 2011) associée à l'utilisation des CRBM</footnote>
Pour ce faire, on a généré entre 50 et 500 échantillons gaussiens. Un niveau de corrélation $\rho$ a été considéré entre deux échantillons successifs. On a alors procédé à l'exécution des tests d'hypothèse cités ci-avant (5000 fois pour chaque cas), et ce, avec un ensemble de 50 et 100 échantillons. Le taux de réjection de l'hypothèse nulle (hypothèse que les échantillons sont bien gaussiens) est reporté sur l'axe des ordonnées ; le niveau de corrélation est reporté en abscisse. On observe la particularité du test KS quant à sa sensibilité au niveau de corrélation. Outre le fait que ce résultat est intéressant à étudier, dans le cas d'une statistique connue, le test KS pourrait être détourné de sa fonction principale, et être utilisé comme détecteur de corrélation. Cette éventualité est laissée en perspective.

<figure><img src="image_3.png" /><figcaption>FIGURE 2.2 – Taux de réjection d'hypothèse de gaussianité d'une distribution à N échantillons corrélés. Le niveau de corrélation (type lag-1) est représenté en abscisse, le taux de réjection est donné par l'axe des ordonnées. Les tests d'hypothèse suivants ont été étudiés : Liliefors, Anderson-Darling (AD), Kolmogorov-Smirnov (KS) ainsi que le $\chi^2$.</figcaption></figure>

## II Régime non surmodé

Le cas du régime surmodé reste une approximation. L'utilisation des chambres s'étend à des régimes où le nombre de modes significatifs n'est pas nécessairement suffisant pour justifier d'une approximation de régime diffus. Un de mes premiers travaux réalisés dans l'hypothèse d'un régime non diffus a porté sur l'étude de la propriété d'ergodicité que l'on retrouve souvent adoptée de manière implicite (en acoustique notamment) sans vraiment de référence permettant de justifier de son domaine de validité. Comme nous le verrons, nos travaux sont basés sur une décomposition du champ de pression en une somme finie d'ondes planes indépendantes. Les autres travaux présentés dans un deuxième temps portant sur l'hypothèse de régime non diffus considèrent quant à eux une décomposition modale.

### II.1 Ergodicité

L'étude des chambres réverbérantes en acoustique est bien antérieure à celles de l'électromagnétisme pour des raisons métrologiques évidentes. Les premières mesures en acoustique ont consisté à relever le temps de réverbération d'une salle ... à l'oreille<fnref n="4" />!

<footnote n="4">4. Sabine WC. The American Architect 1900</footnote>
A notre connaissance, les premiers travaux traitant de l'utilisation du brasseur en acoustique remontent aux années 1950 (BALACHANDRAN 1959), mais son utilisation s'est rapidement faite très rare, au profit d'utilisation de chambres réverbérantes de formes chaotiques dans lesquelles l'approche statistique est spatiale, i.e., une statistique basée sur un échantillonnage spatial du champ. On expose Fig. 2.3 l'intérieur d'une chambre réverbérante acoustique, en l'occurence celle du LNE (Trappes), où l'on peut voir des diffuseurs (plaques gondolés) et des micro suspendus servant à échantillonner le champ de pression. Au centre, une tige rotative munie d'un micro a pour but d'échantillonner de manière dynamique le champ - sachant la taille de la tige, le mouvement rotatif a en fait un effet de brassage sur le champ. En d'autres termes, l'étude de la statistique du champ de pression est spatiale et non sur un ensemble de réalisations (i.e., et non sur un ensemble de positions du brasseur). Il est commun en acoustique de considérer que les statistiques obtenues en prenant un jeu d'échantillons spatiaux et un jeu d'échantillons extraits d'un ensemble de chambres réverbérantes (un échantillon par chambre) convergent : on parle de principe d'ergodicité. Ce principe est bien connu dans l'étude du bruit thermique par exemple où l'on considère que la statistique obtenue via l'ensemble des valeurs d'une mesure unique sur un ensemble de résistances convergera vers la statistique de l'ensemble des valeurs obtenues dans le temps sur une seule résistance. Le but des travaux entrepris a été de fixer un

<figure><img src="image_4.png" /><figcaption>FIGURE 2.3 – Chambre réverbérante acoustique du LNE (Trappes). Les plaques suspendues sont des diffuseurs servant à casser tout degré de symétrie afin d'assurer une chaoticité du milieu ; les microphones suspendus servent quant à eux à échantillonner le champ. Au centre une tige rotative qui sert à échantillonner le champ ; le mouvement de cette tige a un effet de brassage du champ.</figcaption></figure>

cadre de validité de cette hypothèse. Si dans le domaine des signaux aléatoires cette hypothèse est souvent tacite, c'est entre autre parce que le nombre d'échantillons indépendants dans le temps peut être relativement important ; dans le cas des chambres réverbérantes le nombre d'échantillons spatiaux indépendants étant fini, un impact sur la précision des estimateurs peut être limitant. Pour ce faire, le champ de pression (scalaire) considéré, noté *p*, a été exprimé en une somme finie d'ondes planes pour rendre compte de la possibilité de modéliser un régime
non surmodé. Celui-ci a été considéré à une pulsation $\omega$. A notre connaissance, le lien entre le nombre d'ondes planes et le degré de diffusivité du champ n'a pas été établie à l'heure actuelle. De nombreuses questions autour de ce modèle peuvent être intéressantes à approfondir.

Pour mener l'étude, nous avons cherché à analyser la convergence de l'estimateur de moyenne spatiale vers la moyenne d'ensemble vraie. La grandeur utile dans l'étude est l'erreur relative, notée $\epsilon$, et définie telle que (MONSEF et COZZA 2011),

$$ \epsilon = \frac{\hat{\mu}_{p_{\omega}^2} - \mu_{p_{\omega}^2}}{\mu_{p_{\omega}^2}}, \qquad (2.1) $$

où $\hat{\mu}_{p_{\omega}^2}$ est la moyenne spatiale estimée et $\mu_{p_{\omega}^2}$ la moyenne d'ensemble vraie. Afin d'évaluer la précision de l'estimateur spatial, les deux premiers moments de l'erreur ont été dérivés.

Il apparaît que l'estimateur est non biaisé, i.e., la moyenne de l'erreur est nulle. Sa variance en revanche a l'expression suivante (MONSEF et COZZA 2011),

$$ \sigma_{\epsilon}^2 = \left[ \frac{N-1}{N} G_{\Omega}(\mathbf{a}) + \frac{1}{N} \right] (\nu_4 - \nu_2^2), \qquad (2.2) $$

où $\nu_n$ est l'espérance d'ordre $n$ des poids des ondes planes, $G_{\Omega}(\mathbf{a})$ une intégrale faisant inter- venir le lieu géométrique $\mathbf{a}$ inscrit dans un domaine $\Omega$ sur lequel les échantillons spatiaux sont considérés et $N$ correspond au nombre d'ondes planes composant le champ de pression. La densité de probabilité des poids des ondes planes est souvent considérée comme suivant une loi normale. Deux cas sont dès lors intéressants à considérer.

Le premier cas est celui du régime diffus pour lequel $N \mapsto \infty$. La variance suit alors l'évo- lution de $G_{\Omega}(\mathbf{a})$. Il est alors intéressant de constater que suivant la dimensionalité du domaine d'échantillonnage (1D, 2D ou 3D), la vitesse de convergence de l'estimateur diffèrera dans des rapports qui ne suivent pas le rapport des dimensionalités, et ce, en raison de la corrélation spatiale. Pour des échantillons régulièrement espacés, on obtient les décroissances suivantes :

$$ G_{\Omega}(\mathbf{a}) \approx \begin{cases} 0.5/(a/\lambda) & \text{dans le cas 1D} \\ 0.2133/(a/\lambda)^{1.8} & \text{dans le cas 2D} \\ 0.0740(a/\lambda)^{-2} & \text{dans le cas 3D} \end{cases} \qquad (2.3) $$

où $\lambda$ correspond à la longueur d'onde.

Le deuxième cas est celui d'un régime non diffus pour lequel on aurait un échantillonnage tel que $G_{\Omega}(\mathbf{a}) \ll 1/N$. On constate alors que la variance de l'estimateur spatial présenterait une borne basse en $1/N$ liée au degré de mélange du milieu. Ce résultat met en évidence que même si le nombre d'échantillons spatiaux étaient infinis et indépendants, une incertitude plafond existerait au niveau de l'estimateur. On constate que la précision de l'estimateur est conditionnée par le nombre d'ondes planes. Ce résultat est purement théorique et nécessiterait de pouvoir être validé en pratique. La distinction entre les deux termes de la variance de l'estimateur ne s'annonce pas aisée expérimentalement. De plus, pour ce faire, il serait utile d'éclaircir le sens de $N$. Or le nombre $N$ d'ondes planes n'est en réalité pas facile à évaluer à l'heure actuelle et il sera intéressant de voir s'il peut être relié au nombre de degré de liberté du champ. Cette notion de degré de liberté a été étudiée uniquement dans le cadre d'une approche modale (cf. Sec. II.2). Pour un détail des résultats on renvoie à l'article joint en annexe.

Intuitivement le concept d'ergodicité nécessite un niveau de complexité suffisant. Si cette condition est nécessaire, elle n'est pas suffisante. On voit en effet clairement apparaître dans (2.2) le rôle de la corrélation spatiale comme facteur limitant.
## II.2 Approche modale

Dans le domaine de la CEM, la description du champ EM dans une CRBM par une approche modale a été proposée dans un rapport de Lehman (1993) (LEHMAN May 1993). La formalisation du rôle du recouvrement modal et son importance dans la statistique du champ au sein d'une CRBM n'a été publiée qu'en 2011 par Andrea Cozza (COZZA 2011). Pour une meilleure compréhension des explications qui suivront, le champ électrique est exprimée de la manière suivante (COZZA 2011),

$$ \mathbf{E}(\mathbf{r}, f) = \sum_{i=1}^{\infty} \alpha_i \psi_i(f) \mathbf{e}_i(\mathbf{r}), \qquad (2.4) $$

où $\mathbf{e}_i$ est le vecteur propre du ième mode qui n'a aucune forme spécifique supposée; $\alpha_i$ est le poids modal, i.e., la constante de couplage entre le ième mode propre et la source tel que : $\alpha_i = \langle \mathbf{J}_s(\mathbf{r}), \mathbf{e}_i(\mathbf{r}) \rangle_{\Omega_s}$, où $\langle \cdot \rangle_{\Omega_s}$ correspond au produit scalaire appliqué au domaine $(\Omega_s)$ de la source. Ainsi, si la polarisation de la source, i.e., la polarisation du vecteur $\mathbf{J}_s$, est normale au vecteur $\mathbf{e}_i(\mathbf{r})$, le mode ne pourra être excité. La fonction $\psi_i$ est la réponse en fréquence du ième mode. Dans le cadre de milieu à « faibles pertes » cette fonction est approximée à une Lorentzienne.

L'analyse se limitant en un point $\mathbf{r}$, la factorisation modale des modes peut être introduites, telle que,

$$ \mathbf{e}_i(\mathbf{r}) = e_i(\mathbf{r}) \hat{\xi}_i(\mathbf{r}), \qquad (2.5) $$

où $\hat{\xi}_i(\mathbf{r})$ est un vecteur de polarisation supposé uniformément réparti sur $4\pi$ stéradians.

Nous ré-exprimons alors le champ tel que

$$ \mathbf{E}(\mathbf{r}, f) = \sum_{i=1}^{\infty} \tilde{\gamma}_i \psi_i(f) \hat{\xi}_i(\mathbf{r}), \qquad (2.6) $$

où les grandeurs complexes $\tilde{\gamma}_i$ sont considérés comme des poids modaux équivalents, dont la partie réelle et imaginaire suivent respectivement une loi normale. Cette hypothèse n'a aucune justification physique et constitue une étude intéressante en soi. En acoustique, les modèles classiques ne fusionnent pas les poids modaux; les poids modaux côté source et côté récepteur sont considérés comme Gaussiens et indépendants (sous réserve que la source et le récepteur sont suffisamment loin l'un de l'autre). La justification d'une loi normale est donnée par (MCDONALD et KAUFMAN 1988). En électromagnétisme, la grosse différence tient au produit scalaire entre la source et les modes propres. Comme nous allons le voir, le recouvrement entre modes est un des aspects les plus importants dans la statistique du champ; nous expliciterons la grandeur utilisée pour quantifier ce que l'on entend par « statistique du champ ».

### Recouvrement modal

Une illustration du recouvrement entre modes (ou recouvrement modal) présent dans (2.6) est donné Fig. 2.4. Du fait des réponses de type Lorentzienne, on y constate que le nombre de modes qui se recouvrent est en toute rigueur infini!! Ainsi le degré de surmodage ou de diffusivité du champ dans une chambre nécessite intuitivement de définir une bande arbitraire, notée ici $B_e$ sur la Fig. 2.4. Cette bande est centrée sur la fréquence de travail à laquelle la chambre est excitée.

En acoustique des salles, l'effet du recouvrement modal a été approfondi par Schroeder (SCHROEDER et K. KUTTRUFF 1962) qui considérait comme bande $B_e$ la bande à -3 dB
de la réponse d'un mode. En toute rigueur, chaque mode a des pertes différentes car celles-ci découlent directement de la topographie et de la polarisation (dans le cas EM). Toutefois, dans une approche statistique, les modes ne sont, par définition, par explicités de manière déterministe; on en est alors réduit à considérer des pertes moyennes identiques pour des modes « proches ». En électromagnétisme, on quantifie habituellement les pertes par le coefficient de qualité. Celui-ci est alors le facteur de qualité dit *composite* de la chambre; en acoustique il est plutôt d'usage de considérer un coefficient d'amortissement (H. KUTTRUFF 2009) dont l'inverse donne (à un facteur près) une image du temps de réverbération.

La représentation donnée Fig. 2.4 pose alors la question de savoir jusqu'à quel n-ième mode il est nécessaire d'aller pour approcher la statistique du champ observée en pratique. Dit autrement, on s'est intéressé au nombre de bande modale $B_M$ à inclure dans $B_e$. Cette question est traitée en détail dans (MONSEF 2012) jointe en annexe.

On y montre que la réponse n'est pas directe et que la bande $B_e$ à considérer est fonction de la densité modale. Ainsi, à faible recouvrement modal, la bande $B_e$ doit approcher presque 20 fois la bande modale ; ce facteur décroît avec l'augmentation de la densité modale. La justification à ce résultat est que plus le recouvrement modal est faible plus chaque mode aura une influence sur la statistique et il faudra donc aller chercher « loin » chaque mode pouvant avoir une contribution. A contrario, quand la densité modale augmente, le nombre de modes fixant la statistique est d'autant plus vite atteint et permet alors d'avoir un $B_e$ plus étroit.

<figure><img src="image_5.png" /><figcaption>FIGURE 2.4 – (a) Illustration du recouvrement modal (bas à gauche) où chaque mode est représenté. Réponse spectrale résultante (gauche haut). (b) Réponse expérimentale montrant la complexité effective de la superposition modale. $B_e$ est la bande effective à considérer incluant le nombre utile de modes.</figcaption></figure>

### Variabilité de la puissance et statistique du champ

Il est question depuis le début de ce manuscrit de statistique du champ. Il convient de préciser ce que l'on entend par là. Dans le domaine de la CEM et de l'acoustique, la grandeur d'intérêt est le carré de la valeur efficace du champ; en EM, il peut s'agir d'une composante cartésienne du champ. Pour ce dernier cas, la grandeur d'intérêt correspondrait à une image de la puissance
reçue par une (petite) antenne linéairement polarisée. Dans une chambre réverbérante, pour des mesures de puissance rayonnée ou de puissance absorbée, un des objectifs est d'avoir une variabilité la plus réduite possible concernant les grandeurs précitées ci-avant. On se réfère à la variabilité spatiale (sur un tour de brasseur si brasseur il y a) afin de ne pas à avoir à se préoccuper de l'emplacement de l'objet sous test; en d'autres termes, que la mesure puisse être indépendante de la position des objets considérés. On s'intéresse alors à la variabilité de la « puissance ». Ce choix est orienté par les normes de standardisation pour les mesures de puissance rayonnée. On pourrait s'intéresser à la variabilité de la densité d'énergie électrique (ou magnétique) comme dans (COZZA 2011). La notion de « statistique du champ » est donc un abus de langage, notamment en électromagnétisme, la grandeur étant vectorielle. La grandeur considérée est néanmoins intéressante à expliciter puisqu'en régime diffus par exemple, la variabilité de la densité d'énergie est trois fois inférieure à celle de la puissance !

Formellement, si $P_r$ est la puissance reçue, la variabilité, notée $\zeta_P^2$ sera définie telle que,

$$ \zeta_P^2 = \frac{E[P_r^2]}{E^2[P_r]} - 1, \qquad (2.7) $$

où $E[\cdot]$ désigne l'opérateur de moyenne d'ensemble.

On montre que $\zeta_P^2$ est telle que (MONSEF et COZZA 2014),

$$ \zeta_P^2 = 1 + \frac{8}{5\pi M_M}. \qquad (2.8) $$

Ces travaux ont permis d'apporter des résultats à intérêt métrologique car ceux-ci précisent les intervalles de confiance à considérer lorsqu'une mesure est effectuée en CRBM. On montre Fig. 2.5 l'évolution des intervalles de confiance des variabilités estimées $\zeta_P^2$ en fonction du nombre $M_M$ de modes recouverts dans une bande modale $B_M$ à -3 dB.

<figure><img src="image_6.png" /><figcaption>FIGURE 2.5 – Valeur moyenne (trait gras) et limites des intervalles de confiance de $\zeta_P^2$ fonction du nombre moyen $M_M$ de modes recouverts pour différentes valeurs de $N_{sp}$ positions indépendantes de brasseur (reportées sur le graphe). Les limites supérieures et inférieures de l'intervalle de confiance sont au-dessus et en-dessous du trait gras respectivement.</figcaption></figure>

De la même manière, pour une densité d'énergie $W$, nous aurions une variabilité, notée $\zeta_W^2$, telle que (MONSEF, COZZA, MONSEF et al. 2014),

$$ \zeta_W^2 = \frac{1}{3} + \frac{2}{3\pi M_M}. \qquad (2.9) $$
Le cas de l'acoustique a été également abordé en collaboration avec le Laboratoire National d'Essais. On a montré que la variabilité du carré de la valeur efficace du champ de pression $p$ est tel que (MONSEF, COZZA, RODRIGUES et al. 2014) :

$$ \zeta_{p^2}^2 = 1 + \frac{K^2 - 2}{\pi M_M}, \qquad (2.10) $$

où $K$ est le kurtosis des poids modaux; $K = 3$ si les poids modaux sont considérés comme gaussiens.

Ces résultats ont été dérivés sur la base d'une hypothèse d'indépendance entre fréquences propres. Cette hypothèse est importante à souligner car les statistiques d'espacement entre fréquences propres adjacentes a été un domaine de recherche important, notamment depuis l'émergence de la théorie des matrices aléatoires (MEHTA 2004).

### Espacement moyen entre fréquences propres

Pour les cavités non chaotiques (ou intégrables), présentant des symétries, la densité de probabilité de l'espacement entre fréquences propres correspond à une loi de Poisson (LYON 1969) dont l'expression est,

$$ p(x) = \exp(-x) \qquad (2.11) $$

et dont la représentation est rappelée Fig. 2.6 (traits en tirets). En revanche, quand la cavité est chaotique, la densité de probabilité correspond à une loi de Wigner 5 dont l'expression est la suivante,

$$ p(x) = \frac{\pi x}{2} \exp\left(-\frac{\pi x^2}{4}\right) \qquad (2.12) $$

et dont la loi est représentée Fig. 2.6 (trait gras plein).

La différence majeure réside dans l'espacement moyen des modes adjacents. La densité de Poisson met en évidence que dans les cavités non chaotiques la dégénérescence des fréquences propres est très probable, i.e., des fréquences propres peuvent être confondues. Dans le cas de la distribution de Wigner en revanche, on constate que la dégénérescence ne peut avoir lieu. La pente à l'origine témoigne de ce que l'on appelle l'effet de *répulsion* entre fréquences propres. Cette propriété permet une meilleure distribution (spectrale) des réponses en fréquence et réduit les risques de forte de déplétion modale (COZZA 2012). C'est la raison pour laquelle les acousticiens ont opté pour des chambres réverbérantes à forme chaotique. En CEM, bien que les chambres réverbérantes soient, vue de l'extérieur, comme parallépipédiques, l'insertion du brasseur, d'une antenne et d'un objet sous test permet de rompre les symétries et d'avoir un degré de chaoticité au sein de la chambre.

Le travail le plus complet à notre connaissance sur l'effet des distributions de l'espace-ment moyen entre fréquences propres sur la variabilité est celui de Lyon (1969) en acoustique (LYON 1969). Ce travail est digne d'être une référence y compris pour la communauté des CRBM en CEM. L'auteur avait connaissance de la loi de Wigner mais a opté pour la loi dite des plus proches voisins (Nearest Neighbour (NN) en anglais) dont l'expression est la suivante,

$$ p(x) = 2\pi \cdot x \exp(-x\sqrt{2\pi}) \qquad (2.13) $$

<footnote n="5">5. qui correspond à une loi de Rayleigh de paramètre $2/\pi$</footnote>
et dont une représentation est donnée Fig. 2.6 (trait plein fin). Lyon utilise des résultats publiés douze ans avant en physique nucléaire. Ces résultats, proposés par Guverich et Pevsner (GU-REVICH et PEVSNER 1957) proposaient un modèle de statistique mettant en évidence l'effet de répulsion. Lyon cherchait à voir l'effet de la répulsion sur la statistique du champ en chambre réverbérante. Il a donc utilisé cette densité de probabilité ; en utilisant (2.13), on peut parvenir à une expression littérale de la variabilité ; ce qui n'est pas possible avec la loi de Wigner donnée en (2.12).

Bien que la répulsion soit surestimée, l'auteur a su, de manière intéressante, analyser l'effet de cet effet sur les variabilités. Il a montré que la répulsion du modèle *NN* avait un impact (3dB) sur la variabilité par rapport au cas de Poisson à fort recouvrement modal ; en revanche, l'effet de la répulsion est négligeable à faible recouvrement modal [p.560 Ref. (LYON 1969)]. Le modèle de Wigner étant reconnu comme décrivant le plus correctement les systèmes chao-

<figure><img src="image_7.png" /><figcaption>FIGURE 2.6 – Densité de probabilité de Wigner (trait gras plein), de Poisson (tirets), des plus proches voisins (NN) (trait plein fin) et triangulaire (traits en pointillés).</figcaption></figure>

tiques de l'Ensemble Orthogonal Gaussien (système à réversibilité temporelle), une publication (R. WEAVER 1989) de Richard Weaver a eu pour but de combler ce manque. Une analyse détaillée (MONSEF, COZZA, RODRIGUES et al. 2014) de ses travaux nous a amené à remettre en question les résultats proposés, considérant que l'approximation utilisée était trop grossière. Comme expliqué dans la section II de l'article joint en annexe, le modèle proposé revient à un modèle de Poisson avec absence de dégénérescence. Nous avons proposé un modèle où les fréquences propres sont considérées comme indépendantes, uniformément réparties dans des *slots* de largeur $1/m(f_w)$ où $m(f_w)$ est la densité modale à la fréquence $f_w$. Celà revient à considérer une densité de probabilité triangulaire pour l'espacement moyen entre fréquences propres (cf Fig. 2.6 (traits pointillés)). L'intérêt est d'avoir l'effet de répulsion proche de celle de Wigner (pente à l'origine) avec un modèle simple. Pour le calcul des variabilités associées, nous utilisons un principe d'équivalence entre moyenne d'ensemble et moyenne fréquencielle. Pour le détail de la dérivation nous renvoyons respectivement aux références (MONSEF, COZZA, MONSEF et al. 2014) et (MONSEF, COZZA, RODRIGUES et al. 2014) joints en annexe pour l'électromagnétisme et l'acoustique. Il est à noter que le domaine de validité de cette équivalence repose l'hypothèse classique de pertes identiques pour l'ensemble des modes excités à une fréquence donnée.

Il est à noter que notre remise en cause des travaux de Weaver a remis également en cause un article d'un autre auteur J.L-Davy (J. DAVY 1981). Notre publication nous a valu les foudres
de ces deux auteurs dans un commentaire publié en 2015 (J. L. DAVY et R. L. WEAVER 2015).

<figure><img src="image_8.png" /><figcaption>FIGURE 2.7 – Variabilité estimée $\hat{\zeta}_W^2$ (axe y gauche) de la densité d'énergie électrique, fonction du nombre $M_M$ de modes dans le cas où la chambre est vide (a) et chargée (b). Les déviations relatives à la valeur asymptotique de 1/3 sont reportées (axe des ordonnées de droite). Les résultat expérimentaux (symboles gris), ceux de simulation MC (points noirs) et les résultats analytiques sont reportés (trait plein). Les barres d'incertitude verticales traduisent un intervalle de confiance de 95% et les barres horizontales rendent compte de la dispersion que le paramètre $M_M$ autour de sa valeur moyenne.</figcaption></figure>

### Variabilités du champ

La validation de l'approche a été basée sur une comparaison avec des résultats expérimentaux. Pour le cas électromagnétique, la CRBM de 13 m³ (CentraleSupelec) a été utilisée dans deux scénarios : vide et chargée de mousse absorbante. Le deuxième scénario vise à insérer des pertes modifiant ainsi le recouvrement modal. Les résultats expérimentaux (symboles) sont reportés Fig. 2.7. On observe un bon accord avec le résultat analytique rappelé en (2.9). Le modèle modal a également été simulé par simulation Monte Carlo (MC). Le but étant d'une part, de vérifier la consistence des résultats, et d'autre part, de pouvoir évaluer les intervalles de confiance. On constate que ces derniers sont également en bon accord avec les résultats expérimentaux.

Une approche analogue a été entreprise avec une chambre réverbérante acoustique vide et chargée de quelques absorbants. Pour ce faire, nous avons ouvert une collaboration avec le Laboratoire National d'Essais (LNE) basé à Trappes. Les temps de réverbération, classiquement notés $T_{60}$, ont été évalués par bandes de tiers d'octave et les résultats sont reportés Fig. 2.8a dans les cas respectifs d'une chambre vide (cercles) et d'une chambre chargée (carrés).

En utilisant l'approximation de la densité modale donnée par la formule de Weyl, nous avons pu estimer le nombre $M_M$ de modes qui se recouvrent dans une bande modale ; les résultats sont reportés Fig. 2.8b dans les cas respectifs d'une chambre vide (cercles) et chargée (carrés). On y remarque que l'effet des pertes est particulièrement observable à basse fréquence. Comme expliqué précédemment, il est d'usage en acoustique de procéder à un échantillonnage spatial pour estimer la variabilité de la valeur efficace du champ de pression. Toutefois, la méthode d'échantillonnage spatial dans la chambre réverbérante consiste à faire tournoyer un microphone
<figure><img src="image_9.png" /><figcaption>FIGURE 2.8 – (a) Temps de réverbération estimés lorsque la chambre est vide (cercles) et chargée (carrés) (gauche). (b) Estimation du nombre de modes recouverts dans une bande modale dans les cas respectifs où la chambre réverbérante est vide (cercles) et chargée (carrés) (droite).</figcaption></figure>

<caption>FIGURE 2.9</caption>

au bout d'une perche de 1,5 m de long et de diamètre proche du centimètre. Après analyse, cette technique est en fait perturbative au sens d'un brasseur. Pour cette raison nous avons estimer la statistique en utilisant des points fixes dans la chambre. La dispersion a alors été réduite en utilisant un moyennage fréquenciel<fnref n="6" />.

## Nombre effectif de modes

Les recherches se sont poursuivies autour de la compréhension et de l'élaboration de modèles statistiques du champ pour en estimer la variabilité. Ces travaux ont été motivés par le constat suivant. En toute rigueur, quelque soit la fréquence de travail, le nombre de modes excités est mathématiquement infini et rien ne permet de définir clairement la notion de « modes significatifs » sans introduire une bande de fréquence arbitraire ! Le problème est que, suivant la bande choisie, la significativité d'un mode change !

Nous avons donc cherché à définir un modèle dans lequel nous n'aurions pas à relier ce que l'on entend par « mode significatif » au choix arbitraire d'une bande donnée. Pour ce faire, nous avons tout simplement créé un modèle dans lequel nous ne subissons pas les queues des Lorentziennes qui mènent à un recouvrement modal rigoureusement infini.

Le critère choisi est donc de créer un modèle équivalent dans lequel la réponse en fréquence des modes est une porte en fréquence. Pour que le modèle soit équivalent au modèle physique, il doit mener à une variabilité identique à celles présentées dans les parties précédentes. De la satisfaction de ce critère découle la largeur de bande $\pi.B_M$ à assurer pour les réponses « portes ». Nous avons donc à considérer un scénario où le recouvrement modal est comme montré Fig. 2.11. Cette bande est appelée *bande statistique*.

<footnote n="6">6. La plus large bande modale a été estimée autour de 0.5 Hz sur une bande [250 Hz : 1 kHz]. Une décimation par pas de 2 Hz a alors été appliquée.</footnote>
<figure><img src="image_10.png" /><figcaption>FIGURE 2.10 – Variabilités de la valeur quadratique du champ de pression obtenues expérimentalement (trait fin), analytiquement via la relation (2.10) pour $K = 2.7$ et $K = 3$ (trait plein gras et trait en tirets-pointillés respectivement). Les résultats par simulation MC (points noirs) y sont reportés avec les barres d'incertitude associées pour le cas d'une chambre vide (figure gauche) et d'une chambre chargée (figure de droite). Le modèle de Weaver (trait en tirets) y est reporté ainsi que le modèle de Davy/Lyon (traits pointillés) qui suppose pour ce dernier un espace ment entre fréquences propres régi par une loi de Poisson.</figcaption></figure>

Il en découle alors que le nombre significatifs de modes est $\pi$ fois le nombre recouverts dans une bande modale. Une question émerge alors quant à la signification de ces modes.

Si chaque mode en un point de l'espace est une variable aléatoire, exprimer le champ$^7$ comme une somme de modes revient, d'un point de vue statistique, à considérer une somme de variables aléatoires ... pondérées par la réponse en fréquence des modes. Via le modèle équivalent, il apparaît que la pondération est du tout ou rien. Ce modèle équivalent revient donc à voir le champ comme une somme de variables aléatoires (VA) indépendantes et identiquement distribuées (iid). Le nombre $\pi M_M$ de modes fictifs n'est autre que le nombre de VA iid ; en d'autres termes, et par analogie avec une distribution de type $\chi^2$, $\pi M_M$ est le nombre de degré de liberté du champ intervenant dans la variabilité de ce dernier au sein de la chambre réverbérante.

### Approche statistique : critère minimum de pertinence de l'approche

Comme nous l'avons vu, l'approche statistique est assez commode dans les environnements complexes tels que les chambres réverbérantes. L'utilisation idéale d'une chambre est en régime diffus ou surmodé sachant que ce régime permet une dispersion minimale quant aux performances moyennes attendues.

Dans les normes existantes, il est néanmoins possible et prévu d'utiliser les chambres dans des régimes non diffus. Dans ces cas là, de plus grandes barres d'incertitude sont attendues. Une question est alors de savoir jusqu'à quelle fréquence minimale est-on susceptible d'utiliser une chambre. La fréquence est alors connue sous l'acronyme de LUF pour *Lowest Useable Frequency* (Fréquence d'utilisation Minimale). La définition actuelle de la LUF n'est pas basée sur un critère en lien avec la physique de la chambre.

Sans rentrer dans le débat de la définition de la LUF, le concept de degré de liberté, introduit plus haut, nous permet de proposer une borne basse et un critère pour lesquelles l'approche

<footnote n="7">7. en fait la puissance, la densité d'énergie ou le carré de la valeur rms du champ de pression;-)</footnote>
<figure><img src="image_11.png" /><figcaption>FIGURE 2.11 – Illustration du recouvrement modal dans le cas où les réponses en fréquences sont des Lorentziennes (traits gras) et dans le cas où les réponses en fréquence sont des portes en fréquence (traits fins). Si la bande modale à -3dB est notée $B_M$, on montre que une bande $\pi B_M$ des fonction portes mène à une variabilité identique.</figcaption></figure>

statistique peut faire sens. Comme présenté dans (MONSEF et COZZA 2015a) donné en annexe, le critère qui nous a paru le plus intuitif est d'avoir un seul degré de liberté au minimum, et ce, afin de pouvoir considérer au moins un terme dans la décomposition modale liée au modèle équivalent introduit plus haut. L'avantage de cette borne est son adaptation au degré de pertes présent dans la chambre. Ainsi, dans le cas d'une chambre avec un objet à faible perte, ou au contraire, une chambre chargée par un objet dissipatif, la fréquence associée à un degré unique de liberté différera. Un pire cas pourra être associé au cas d'une chambre vide. Le lien avec la LUF dépendra des barres d'incertitude acceptée pour la variabilité. Ceci n'est pas explicité dans les normes actuelles (Reverberation Chamber Test Methods 2011). Dans le contexte des chambres, on rappelle que c'est la statistique du champ maximum qui intéresse la communauté de la CEM; l'étude du maximum du champ avec une approche modale n'est pas aisée et n'a pas été développée à notre connaissance. En revanche, un des aspects intéressant la communauté est de pouvoir proposer un modèle qui puisse décrire le cas d'une chambre réverbérante du cas sous-modé jusqu'au régime asymptotique.

## III Du Régime sous-modé au régime sur-modé

Une des questions pendantes est de pouvoir proposer un modèle statistique d'une des grandes d'intérêt et de le relier à la physique de la chambre; la contrainte étant de pouvoir décrire n'importe quel régime de fonctionnement avec un tel modèle. Une approche statistique basé sur la distribution de Weibull, sans lien avec la physique de la chambre, a été proposé par Bunting (BUNTING 2002); le but étant de pouvoir utiliser une loi de distribution de la puissance mesurée par une antenne pour n'importe quel régime de fonctionnement de la chambre. Cette distribution, notée ici $f_X(x)$, a l'expression suivante,

$$f_X(x) = ab^{x-1} \exp(-ax^b), \qquad (2.14)$$
<figure><img src="image_12.png" /><figcaption>FIGURE 2.12 – Estimation du nombre $M_M$ de modes se recouvrant dans une bande modale en fonction de la fréquence pour une chambre réverbérante vide de 13m³; le modèle de Weyl y est reporté également (trait plein)</figcaption></figure>

où $a$ et $b$ sont respectivement les paramètres d'échelle et de forme de la distribution. Dans le cas d'un régime surmodé, le champ étant Gaussien, la puissance suit une loi exponentielle, i.e., $a = b = 1$.

Comme on l'a vu précédemment, un des paramètres physiques, pour ne pas dire le seul, mesurant le degré de « surmodage » de la chambre est le nombre de modes significatifs découlant du recouvrement modal; sans perte de généralité, le nombre $M_M$ de modes se recouvrant dans la bande modale moyenne peut également être considérée sachant un facteur $\pi$ constant entre les deux grandeurs. Notre proposition a donc été de relier l'approche modale à la distribution de Weibull afin d'exprimer son paramètre de forme en fonction du nombre $M_M$ de modes.

L'approche proposée a consisté en une identification des variabilités obtenues par le modèle de Weibull et l'approche modale. Comme montré dans (MONSEF et COZZA 2015b), celà revient à étudier l'égalité suivante,

$$ \frac{\mathrm{Var}[P]}{\mathrm{E}^2[P]} = 1 + \frac{8}{5\pi M_M} = \frac{\Gamma\left(1 + \frac{2}{b}\right)}{\Gamma^2\left(1 + \frac{1}{b}\right)} - 1, \qquad (2.15) $$

où $\Gamma(\cdot)$ est la fonction d'Euler.

Comme on peut le constater, une expression reliant le paramètre de forme $b$ au nombre $M_M$ de modes ne peut être explicitée. Par résolution numérique, on obtient l'évolution donnée Fig. 2.13; on propose dans la figure une courbe de tendance reliant le paramètre statistique au paramètre physique.

Pour valider cette relation, on a procédé à des estimations du paramètre $b$ par le biais de relevés expérimentaux. Ces mesures ont consisté à faire des mesures sur une bande [700 MHz :3 GHz] avec 1000 points en fréquence pour chaque composante cartésienne du champ; l'opération a été effectuée pour 100 positions de brasseur. On a alors regroupé les échantillons par 5 dans l'espace des fréquences. L'estimation a ainsi été effectuée sur la base de 1500 échantillons pour les 200 fréquences considérées.

Par l'estimation de la variabilité, on a pu remonter au nombre $M_M$ de modes qui se recouvrent dans une bande modale. On reporte les résultats Fig. 2.12. On peut constater que la chambre est sous-modée au sens de Schroëder pour des fréquences inférieures à 1.5 GHz approximativement. On a établi une courbe de tendance (trait plein) afin de pouvoir établir dans un second temps une comparaison avec les valeurs estimées de $b$.
<figure><img src="image_13.png" /><figcaption>FIGURE 2.13 – Evolution du paramètre de forme b de la distribution de Weibull en fonction du nombre $M_M$ de modes se recouvrant. On propose également une expression de la courbe de tendance afin d'avoir une relation explicite entre les deux paramètres.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>FIGURE 2.14 – Evolution du paramètre de forme b estimé en fonction de la fréquence (gauche) et du nombre $M_M$ (droite). L'estimation s'est faite sur la base de 1500 échantillons. Les symboles gris ont été obtenus expérimentalement, la courbe en trait plein découle de la relation obtenue Fig. 2.13 couplée à l'estimation du nombre $M_M$ en fonction de la fréquence. Les barres d'incertitude ont été obtenues par simulation Monte Carlo.</figcaption></figure>

A l'aide des données expérimentales, nous avons également estimé les fonctions de distribution de la puissance, et par ce biais, nous avons pu estimer les facteurs de forme b de la loi de Weibull. On peut voir les ordres de grandeur à l'aide de la Fig. 2.14. Comme attendu une convergence vers la valeur unité est observée avec l'augmentation de la fréquence. A l'aide de la courbe de tendance $M_M(f)$ obtenue via la Fig. 2.12, nous avons déduit la courbe en trait plein de la Fig. 2.14 gauche. Une comparaison avec l'évolution de $b(M_M)$ est également présentée Fig. 2.14 (droite). Afin de savoir si la dispersion des points expérimentaux (cercles gris) était à imputer à l'incertitude statistique, nous avons procédé à une simulation du modèle modal en un point. Nous obtenons des barres d'incertitude en bon accord avec la dispersion observée expérimentalement, confirmant que cette dispersion est liée à la nature aléatoire du champ EM au sein de la chambre.
# IV Conclusion et perspectives

En conclusion, la plupart des travaux ont porté sur la statistique de la moyenne du champ au carré au sein d'une chambre réverbérante. Les travaux menés sur cette thématique ont permis de mieux estimer les intervalles de confiance en fonction du degré de surmodage d'une chambre. Une des contributions jugée comme intéressante concerne la possibilité d'associer un nombre de degrés de liberté suivant la statistique du champ au sein de l'enceinte.

Une des perspectives, déjà entamée par ailleurs, porte sur l'extension et l'exploitation de ces résultats au domaine des télécommunications. Les chambres réverbérantes étant vues comme des émulateurs de canaux, les résultats obtenus permettent d'étendre l'utilisation de ces dernières à de nouveaux types de *fading*.

Toutefois, un des défis sur la thématique des chambres réverbérantes porte sur la statistique des maximas de champ. Dans les tests de susceptibilité électromagnétique ce sont les maxima qui intéressent la communauté et dans des conditions de surmodage intermédiaire, la statistique de ces maxima n'est pas aisée à trouver.

Cette étude est laissée en perspective. Pour y parvenir, une des étapes porterait sur l'étude de la statistique des poids modaux. A défaut de pouvoir statuer sur ce point, la gaussianité de ces poids est considérée jusque lors mais rien n'indique que ce modèle soit le bon. Un travail en ce sens a été entamé. Avec le modèle modal équivalent proposé il serait intéressant de voir si la statistique des maxima serait alors vérifiée expérimentalement.

Un autre point laissé en perspective concerne la décomposition du champ en ondes planes avec notamment l'étude d'une corrélation potentielle entre un nombre fini d'ondes planes qu'inclurait le modèle et le recouvrement modal.

# Références

BALACHANDRAN, C. (1959). « Random sound field in reverberation chambers ». *J. Acoust. Soc. Am* 31, p. 1319 (cf. p. 9).

BUNTING, C. F. (2002). « Statistical characterization and the simulation of a reverberation chamber using finite-element techniques ». *IEEE Trans. Electromagn. Compat.* 44.1, p. 214-221 (cf. p. 19).

COOK, R. K., R. WATERHOUSE, R. BERENDT, S. EDELMAN et M. THOMPSON JR (1955). « Measurement of correlation coefficients in reverberant sound fields ». *J. Acoust. Soc. Am.* 27.6, p. 1072-1077 (cf. p. 7).

COZZA, A. (2012). « Probability Distributions of Local Modal-Density Fluctuations in an Electromagnetic Cavity ». *IEEE Trans. Electromagn. Compat.* 54.5, p. 954-967 (cf. p. 7, 14).

COZZA, A. (2011). « The role of losses in the definition of the overmoded condition for reverberation chambers and their statistics ». *IEEE Trans. Electromagn. Compat.* 53.2, p. 296-307 (cf. p. 11, 13).

DAVY, J. (1981). « The relative variance of the transmission function of a reverberation room ». *J. Sound Vib.* 77.4, p. 455-479. ISSN: 0022-460X (cf. p. 15).

DAVY, J. L. et R. L. WEAVER (2015). « Comment on "Relative variance of the mean squared pressure in multimode media : Rehabilitating former approaches" | *J. Acoust. Soc. Am.* 136, 2621-2629 (2014)». *J. Acoust. Soc. Am.* 137.3, p. 1598-1601 (cf. p. 16).

GUREVICH, I. et M. PEVSNER (1957). « Repulsion of nuclear levels ». *Nuclear Physics* 2.5, p. 575-581 (cf. p. 15).
HILL, D. (1998). « Plane-wave integral representation for fields in reverberation chambers ». *IEEE Trans. Electromagn. Compat.* 40, p. 209-217 (cf. p. 6).

KUTTRUFF, H. (2009). « Room acoustics ». 5th. Oxon : Spon Press. Chap. 3 (cf. p. 12).

LEHMAN, T. H. (May 1993). « A Statistical Theory of Electromagnetic Fields in Complex Cavities ». *Interaction Notes, USAF Phillips Laboratory* Note 494 (cf. p. 11).

LYON, R. (1969). « Statistical analysis of power injection and response in structures and rooms ». *J. Acoust. Soc. Am* 45, p. 545 (cf. p. 14, 15).

MCDONALD, S. et A. KAUFMAN (1988). « Wave chaos in the stadium : statistical properties of short-wave solutions of the Helmholtz equation ». *Phys. Rev. A* 37.8, p. 3067 (cf. p. 11).

MEHTA, M. L. (2004). *Random Matrices*. Sous la dir. d'E. INC. Third. T. 142. Amsterdam, Oxford, San Diego, London : Academic Press (cf. p. 14).

MONSEF, F. (2012). « Why a Reverberation Chamber Works at Low Modal Overlap ». *IEEE Trans. Electromagn. Compat.* 54.6, p. 1314-1317 (cf. p. 12).

MONSEF, F., A. COZZA, F. MONSEF et A. COZZA (2014). « Average number of significant modes excited in a mode-stirred reverberation chamber ». *IEEE Trans. Electromagn. Compat.* 56.2, p. 259-265 (cf. p. 13, 15).

MONSEF, F. et A. COZZA (2011). « Limitations of the equivalence between spatial and ensemble estimators in the case of a single-tone excitation ». *J. Acoust. Soc. Am.* 130, p. 1943 (cf. p. 10).

MONSEF, F. et A. COZZA (2014). « Variability and confidence intervals of the power measured in a reverberation chamber ». *IEEE Trans. Electromagn. Compat.* 56.5, p. 1238-1241 (cf. p. 13).

MONSEF, F. et A. COZZA (2015a). « A Possible Minimum Relevance Requirement for a Statistical Approach in a Reverberation Chamber ». *IEEE Trans. Electromagn. Compat.* 57.6, p. 1728-1731 (cf. p. 19).

MONSEF, F. et A. COZZA (2015b). « Predicting Power Distribution Inside a Reverberation Chamber ». *URSI-RASC Conference Proceedings*. IEEE (cf. p. 20).

MONSEF, F., A. COZZA, D. RODRIGUES, P. CELLARD et J.-N. DUROCHER (2014). « Relative variance of the mean-squared pressure in multimode media : Rehabilitating former approaches ». *J. Acoust. Soc. Am.* 136.5, p. 2621-2629 (cf. p. 14, 15).

MONSEF, F., R. SERRA et A. COZZA (2015). « Goodness-of-Fit Tests in Reverberation Chambers : Is Sample Independence Necessary? » *IEEE Trans. Electromagn. Compat.* 57.6, p. 1748-1751 (cf. p. 7).

MORSE, P. M. et R. H. BOLT (1944). « Sound waves in rooms ». *Rev. Mod. Phys.* 16.2, p. 69 (cf. p. 6).

*Reverberation Chamber Test Methods* (2011). International Electrotechnical Commission (IEC), Std. 61 000-4-21 (cf. p. 7, 19).

SCHROEDER, M. R. et K. KUTTRUFF (1962). « On frequency response curves in rooms. Comparison of experimental, theoretical, and Monte Carlo results for the average frequency spacing between maxima ». *J. Acoust. Soc. Am.* 34.1, p. 76-80 (cf. p. 7, 11).

WEAVER, R. (1989). « On the ensemble variance of reverberation room transmission functions, the effect of spectral rigidity ». *J. Sound Vib.* 130.3, p. 487-491. ISSN : 0022-460X (cf. p. 15).

WEYL, H. (1912). « Das asymptotische Verteilungsgesetz der Eigenwerte linearer partieller Differentialgleichungen (mit einer Anwendung auf die Theorie der Hohlraumstrahlung)(Asymptotic distribution of the partial differential equations eigenvalues (with an application to black body theory)) ». *Mathematische Annalen* 71.4, p. 441-479. ISSN : 0025-5831 (cf. p. 6).
# Limitations of the equivalence between spatial and ensemble estimators in the case of a single-tone excitation

Florian Monsef<sup>a)</sup> and Andrea Cozza

Département de Recherche en Électromagnétisme, L2S, UMR8506

Univ Paris-Sud, SUPELEC, CNRS

3 rue Joliot Curie - 91190 Gif-sur-Yvette, France

(Dated: August 14, 2011)

The ensemble-average value of the mean-square pressure is often assessed by using the spatial-average technique, underlying an equivalence principle between spatial and ensemble estimators. Using the ideal-diffuse-field model, the accuracy of the spatial-average method has been studied theoretically forty years ago in the case of a single-tone excitation. This study is revisited in the present work on the basis of a more realistic description of the sound field accounting for a finite number of plane waves. The analysis of the spatial-average estimator is based on the study of its convergence rate. Using experimental data from practical examples, it is shown that the classical expression underestimates the estimator uncertainty even for frequencies greater than Schroeder's frequency, and that the number of plane waves may act as lower bound on the spatial-average estimator accuracy. The comparison of the convergence rate with an ensemble-estimator shows that the two statistics cannot be regarded as equivalent in a general case.

PACS numbers: 43.55.Cs, 43.55.Br, 43.55.Gx

## I. INTRODUCTION

In room acoustics, the interest of statistical approaches lies in the ability of predicting statistical properties of the sound field independently of the fine details of a specific configuration. Indeed, the configuration of a room and the position of the source impact on the modal topographies or equivalently on the directions and the amplitudes of the plane waves composing the sound field. If one considers an ensemble of rooms with approximately the same volume and the same reverberation time, one should recover the same estimated moments in accordance with the statistical model.

Statistical models have been developed for more than a half century. The most common one, called the ideal-diffuse-field model, is based on random wave theory¹. The model considers, on the one hand, the reverberant part of the sound field, that is the sound field measured at some distance from the source and, on the other hand, that the modal overlap is sufficiently great. It has been shown that a minimum average overlap of three modes is needed², from which ensues Schroeder's frequency. Other studies³⁻⁵, based on a modal approach, aimed at deriving models that could be applied to a more general class of systems, i.e., to configurations where a low modal overlap can be assumed.

In practice, the validation of any model depends on estimation techniques based on the use of samples from which estimates of the statistical moments (generally the first two moments) must be derived. The random samples can be obtained mainly by two different methods. The first method consists in sampling the sound pres-

sure in the central volume of a room<sup>4,6,7</sup>, for a given configuration of the sound-field topography. The estimation is based on a single realization of the sound-field topography. The estimation method using spatial samples will hereafter be referred to as the spatial-estimator method. The second method consists in using a single spatial sample for an ensemble of different sound-field topographies<sup>8,9</sup> that can be obtained electronically, by the use of a varying tone signal<sup>10</sup>, or mechanically, using a large rotating vane<sup>8,9</sup>. This latter is particularly well-suited for single-tone excitation which is the case of interest in the present study. Hereafter, this estimation method will be referred to as an ensemble-estimator method.

Using the squared value of the sound pressure related to each samples set, two estimations of the mean-square pressure can then be computed. Although never demonstrated, to the best of our knowledge, the estimates resulting from the spatial- and ensemble-average estimators are often regarded as equivalent. This equivalence principle implies a similar asymptotic convergence rate of both estimators. The equivalence assumption is directly related to the mixing (and consequently ergodic) property of dynamical systems<sup>11</sup>. Using a classical wave theory approach, an analogous mixing property is found within the ideal-diffuse-field description<sup>12</sup>. This equivalence assumption has been widely used to optimize spatial-averaging techniques, studied in detail by Schroeder<sup>13,14</sup>, Lubman<sup>15-17</sup> and Waterhouse<sup>18</sup> forty years ago, to provide the effective number of discrete independent spatial samples for practical experiments. However, in the case of a single-tone excitation, these works are restricted to frequencies above Schroeder's frequency.

The aim of the present paper is to study theoretically, in the case of a harmonic excitation and considering only the reverberant part of the sound pressure, whether the spatial- and ensemble-average estimators can be regarded

<footnote><sup>a)</sup>florian.monsef@lss.supelec.fr</footnote>
as equivalent. To this end, considering the mean-square pressure, expressions of the first two moments of the estimators are derived. The ensemble and spatial estimators are shown to be both unbiased but with different convergence rates. As a result, the confidence intervals related to both estimators cannot be regarded as identical.

To derive these results, a plane-wave model ensuing from a modal expansion is used to express the sound pressure. Recalling on the one hand that a mode can be expanded into a strictly finite number of plane waves¹⁹, and considering on the other hand a finite number of overlapping modes, the model used includes a discrete number of plane waves. Hence, this approach ensures a realistic description of the field and is presented in section II. Using this sound-field model, we establish in section III a general expression of the spatial-average estimator variance. The accuracy of the estimation is assessed on the basis of its convergence rate, and it is shown to present a lower bound directly related to the number of plane waves composing the sound pressure. In section IV we relate the number of plane waves to the room parameters, providing a predictive model applied to practical-room examples. We show that the lower bound of the convergence rate has a significant impact on the accuracy of the estimated value even for frequencies above Schroeder's frequency.

## II. HARMONIC SOUND FIELD MODEL

### A. From modal theory to discrete plane-wave expansion

In any enclosure with reflective boundaries (room, theater, reverberation chamber,...), a sound pressure induced by a source signal is made up of standing waves related to the normal modes that form a complete set of orthogonal functions¹⁹. Accordingly, a harmonic sound field can be expressed as a sum of M dominant overlapping modes at a given frequency. Considering a harmonic excitation, the sound pressure will be expressed in the frequency domain. The resulting complex pressure $p_\omega(\mathbf{r})$ can then be stated as

$$p_{\omega}(\mathbf{r}) = \sum_{m=1}^{M} \alpha_{m} p_{m}(\mathbf{r}) \psi_{m}(f), \quad (1)$$

where $p_m(\mathbf{r})$, $\alpha_m$ and $\psi_m(f)$ correspond to the modal topography, the modal weight and the frequency response of the $m^{th}$-mode, respectively.

As it will be presented thereafter, the spatial-average estimator consists basically in computing a spatial integral in which the modal topographies are needed. The computation is achievable if explicit expressions of the modes are available which is only the case in very specific configurations. This points out a serious limitation of the modal expansion in (1) as a way to study the accuracy of the spatial-average estimator.

Nevertheless, in the case of a standard room, a mode can be expanded into a finite number of plane waves. Consequently, for a finite modal overlap, the sound field can be rigourously expressed by means of a finite number

of plane waves. To the best of our knowledge, this is the first analysis in room acoustics introducing such a model, allowing a more realistic description of the sound field. Indeed, it has been known for a long time that the infinite number of plane waves assumed in the ideal diffuse model is an unrealistic concept²⁰ and is not physically justified. In practical terms, the finite number of plane waves accounts for the finite number of possible paths (as illustrated in Ref.10), or periodic orbits²¹, inducing cavity resonances and stationary waves. The resulting complex pressure $p_\omega(\mathbf{r})$ expressed in (1) can then be restated as

$$p_{\omega}(\mathbf{r}) = \sum_{p=1}^{N} \gamma_{p} e^{-jk_{0} \hat{k}_{p} \cdot \mathbf{r}}, \quad (2)$$

where N is the number of plane waves, $\mathbf{r}$ is the position at which the pressure is sampled and $k_0 = \omega/c$ is the propagation constant for the homogeneous medium filling the cavity, with $\omega$ the angular frequency of the harmonic excitation signal and c the sound speed. The $p^{th}$-plane wave propagates along the direction $\hat{k}_p$, with a complex amplitude $\gamma_p$.

We want to stress that the plane-wave model introduced in (2) should be regarded as a modal approach expressed in the Fourier basis. In other words, using (1) or (2) is totally equivalent. The main advantage of (2) over (1) however, is that the exponential term provides an explicit spatial dependence well suited to the study in section III of the spatial estimator. This will allow thereafter a rigorous derivation of the statistical moments of the spatial estimator.

### B. Ensemble statistics

In practice, the plane-wave sets of parameters $\{\hat{k}\} = \{\hat{k}_1, \dots, \hat{k}_N\}$ and $\{\gamma\} = \{\gamma_1, \dots, \gamma_N\}$ cannot be measured, the only observable being the pressure. As commonly done in statistical mechanics and in the modeling of complex systems, the use of (2) allows studying sound field statistics in a simple way, as soon as the two plane-wave sets of parameters $\{\hat{k}\}$ and $\{\gamma\}$ are treated as random variables.

Hence, to make the derivation of our approach simpler, we adopt the following assumptions on $\{\hat{k}\}$ and $\{\gamma\}$

1. the two sets will be assumed to be independent

2. for each set, variables identified by different indexes $p$ in (2) will be regarded as independent and identically distributed (iid).

Although the expression of the sound pressure in (2) holds in a general case, i.e., for the sound field near the source and at some distance from the source, the hypotheses made on the sets $\{\hat{k}\}$ and $\{\gamma\}$ implies physically a certain degree of diffusion of the sound field. Accordingly, the model that will be established in the present study is meaningful only for the reverberant part of the field with a certain degree of diffusion.
The mean-square pressure being the quantity of interest, we use (2) to express the square of the pressure as follows

$$p_{\omega}^{2}(\mathbf{r}) = \sum_{p=1}^{N} |\gamma_{p}|^{2} + \sum_{i \neq j}^{N} \gamma_{i} \gamma_{j}^{*} e^{-jk_{0}(\hat{k}_{i} - \hat{k}_{j}) \cdot \mathbf{r}}, \quad (3)$$

where * stands for the complex conjugate.

In order to derive the statistical properties of (3), we start by expressing its ensemble-average value $\mu_{p_{\omega}^2}(\mathbf{r})$ defined as

$$\mu_{p_{\omega}^2}(\mathbf{r}) = E[p_{\omega}^2(\mathbf{r})] = \sum_{p=1}^{N} E[|\gamma_p|^2] + \sum_{i \neq j}^{N} E[\gamma_i \gamma_j^*] E[e^{-jk_0(\hat{k}_i - \hat{k}_j) \cdot \mathbf{r}}], \quad (4)$$

where $E[\cdot]$ stands for the ensemble-average operator. It is noteworthy that we have taken advantage of the independence assumption between the sets $\{\gamma\}$ and $\{\hat{k}\}$ to split the ensemble operator in the last term in (4).

One of the properties sought in room acoustics applications, e.g., power measurements, is the statistical uniformity of the sound field. In other words, $\mu_{p_{\omega}^2}$ is expected to be independent of the position $\mathbf{r}$. One way to justify this property in some previous work$^{14}$, is to consider $\{\gamma\}$ as a set of zero mean complex random variables, such as

$$E[\gamma_i^* \gamma_j] = E[\gamma_i^*] E[\gamma_j] = 0 \quad \forall i \neq j. \quad (5)$$

In fact, this assumption is not necessary, as long as the sampling position is in a region where the ensemble of the plane waves cover uniformly angles of arrival over $4\pi$ steradian. It is worth noting that this property holds even for a finite number of plane waves. Assuming the $\{\hat{k}\}$ to be iid leads to have the ensemble-average value of the complex exponential term in (4) to be equal to zero. Therefore,

$$\mu_{p_{\omega}^2}(\mathbf{r}) = \sum_{p=1}^{N} E[|\gamma_p|^2] = N E[|\gamma_p|^2]. \quad (6)$$

### III. SPATIAL ESTIMATOR

In order to characterize in detail the properties of this estimation method, we make use, for convenience, of a continuous spatial average which provides the following estimated value $\hat{\mu}_{p_{\omega}^2}$

$$\hat{\mu}_{p_{\omega}^2} = \frac{1}{M_{\Omega}} \int_{\Omega} p_{\omega}^2(\mathbf{r}) \, d\mathbf{r} = \langle p_{\omega}^2(\mathbf{r}) \rangle_{\Omega}, \quad (7)$$

where $\langle \cdot \rangle_{\Omega}$ stands for the space averaging operator as taken over the sample region $\Omega$.

In order to derive a more explicit expression of the spatial-average estimation, we substitute (3) into (7).

Accordingly,

$$\begin{aligned} \hat{\mu}_{p_{\omega}^2} &= \frac{1}{M_{\Omega}} \sum_{i=1}^{N} \sum_{j=1}^{N} \gamma_i \gamma_j^* \int_{\Omega} e^{-jk_0(\hat{k}_i - \hat{k}_j) \cdot \mathbf{r}} \, d\mathbf{r} \\ &= \sum_{i=1}^{N} \sum_{j=1}^{N} \gamma_i \gamma_j^* \langle e^{-jk_0(\hat{k}_i - \hat{k}_j) \cdot \mathbf{r}} \rangle_{\Omega}. \end{aligned} \quad (8)$$

Organizing the set of the plane-wave weights $\{\gamma\}$ into the $\gamma$ column vector yields the following quadratic form

$$\hat{\mu}_{p_{\omega}^2} = \gamma^H \mathbf{P} \gamma, \quad (9)$$

where $\mathbf{H}$ stands for the Hermitian operator and $\mathbf{P}$ is a matrix with elements $P_{ij}$ defined as

$$P_{ij} = \langle e^{-jk_0(\hat{k}_i - \hat{k}_j) \cdot \mathbf{r}} \rangle_{\Omega}. \quad (10)$$

In the case of a 3D cartesian sampling volume, the elements of $\mathbf{P}$ can be expressed as,

$$P_{ij} = \frac{1}{M_{\Omega}} \int_{\Omega} e^{-jk_0(\hat{k}_i - \hat{k}_j) \cdot (x\hat{n}_1 + y\hat{n}_2 + z\hat{n}_3)} \, dx \, dy \, dz, \quad (11)$$

where $\hat{n}_1$, $\hat{n}_2$ and $\hat{n}_3$ are the Cartesian unit vectors related to the $x$, $y$ and $z$ components, respectively. If we consider a square parallelepiped whose sides are referred to as $a_1$, $a_2$, $a_3$, respectively, and if the origin is taken at the middle of the sampling volume $\Omega$, we can easily restate (11) as

$$P_{ij} = \prod_{v=1}^{3} \mathrm{sinc} \left( \pi \frac{a_v}{\lambda} (\hat{k}_i - \hat{k}_j) \cdot \hat{n}_v \right), \quad (12)$$

where sinc(·) is the sine cardinal function, $\lambda$ is the wavelength in the filling medium. We note that the elements of matrix $\mathbf{P}$ are a function of the random variables $\hat{k}_i - \hat{k}_j$.

Beyond the mathematical convenience provided by the concise form of (9), the matrix $\mathbf{P}$ is of particular interest in understanding the influence of the non-diagonal terms of $\mathbf{P}$ on the estimator accuracy. To illustrate our statement, suppose that $\mathbf{P}$ reduces to the identity matrix $\mathbf{I}$. A perfect estimator is obtained since $\hat{\mu}_{p_{\omega}^2}$ converges to the ensemble average given in (6). The scenario ensuring this convergence is however unrealistic since corresponding to an infinite sampling volume as illustrated by relation (12). Indeed, if $a_v$ tends to infinity the cross-terms $P_{ij}$ tend to zero.

In practice however, the sampling volume being finite, non-diagonal terms of $\mathbf{P}$ will appear, i.e., $\mathbf{P} \neq \mathbf{I}$, inducing inevitably estimation errors. In other words, the analysis of the spatial-estimator accuracy lays on the analysis of the eigenvalues of $\mathbf{P}$.

In order to study the accuracy of the estimate $\hat{\mu}_{p_{\omega}^2}$, we need to compute the moments of the spatial estimator with respect to ensemble statistics.

### A. Statistics of the error induced by a spatial estimator
The accuracy of the spatial estimator lies in the ability of the corresponding estimate to provide the true ensemble-average value. Accordingly, the quantity

$$ \epsilon = \frac{\hat{\mu}_{p_{\omega}^2} - \mu_{p_{\omega}^2}}{\mu_{p_{\omega}^2}} \qquad (13) $$

provides a measure of the relative error of the spatial-average value with respect to the ensemble-average value. In order to assess the accuracy of the spatial estimator the two statistical moments of (13) must be studied. The first moment, i.e., the mean value, informs on the asymptotic convergence of the spatial estimator towards the ensemble-mean value, i.e., informs on the biased or unbiased nature of the estimator. The second moment, i.e., the variance, is crucial since it assesses the confidence interval of the estimated value, i.e., indicates the accuracy of the estimate. Before proceeding, we express the relative error $\epsilon$ in a different form for a more convenient mathematical derivation. We start by expressing (13) in a more explicit form. Using (6) and (9) yields

$$ \epsilon = \frac{\hat{\mu}_{p_{\omega}^2} - 1}{\mu_{p_{\omega}^2}} = \frac{\gamma^H P \gamma}{N E [|\gamma_p|^2]} - 1. \qquad (14) $$

It is convenient to introduce the random vector

$$ \mathbf{v} = \frac{\gamma}{\sqrt{E[|\gamma_p|^2]}}, \qquad (15) $$

which is made up of iid normalized random variables, i.e., $E[|v_i|^2] = 1, \forall i \in [1, N]$. That allows recasting (14) in a simpler form

$$ \epsilon = \frac{1}{N} \mathbf{v}^H \mathbf{P} \mathbf{v} - 1. \qquad (16) $$

As explained in the previous subsection, a convenient way of analyzing the accuracy of the spatial estimator is to introduce the eigenvalues of $\mathbf{P}$ using the following change of basis

$$ \mathbf{P} = \mathbf{X} \mathbf{\Lambda} \mathbf{X}^{-1}, \qquad (17) $$

where $\mathbf{\Lambda}$ is a diagonal matrix composed of the eigenvalues $\lambda_i$ of $\mathbf{P}$. Consequently,

$$ \epsilon = \frac{1}{N} \mathbf{w}^H \mathbf{\Lambda} \mathbf{w} - 1, \qquad (18) $$

where $\mathbf{w} = \mathbf{X}^{-1}\mathbf{v}$ is the vector of the plane-wave weights in the new basis. $\mathbf{X}$ being a unitary matrix, the probability density function of the random variables $v_i$ and $w_i$ is the same. Moreover, by computing the covariance matrix of the set $\{w\}$, it can easily be shown that the statistical properties assumed for the $v_i$ random variables are conserved for the $w_i$.

## B. Mean of the error

The statistical description of the relative error deals with a multivariate statistics analysis. To perform such

analysis we need to consider a joint probability density function $p(\{v\}, \{\hat{k}\})$. Recalling the independence assumptions between the two sets of random parameters $\{v\}$ and $\{\hat{k}\}$, the joint probability function can be factored into the following way

$$ p(\{v\}, \{\hat{k}\}) = p_v(\{v\}) p_{\hat{k}}(\{\hat{k}\}), \qquad (19) $$

where $p_v(\{v\})$ and $p_{\hat{k}}(\{\hat{k}\})$ are the marginal probability density function of the random sets $\{v\}$ and $\{\hat{k}\}$, respectively.

Hereafter, we will adopt the notation $E_v[\epsilon]$ related to the conditional average of the random function $\epsilon$ conditioned by the random set $\{v\}$, defined as

$$ E_v[\epsilon] = \int \epsilon (\{v\}, \{\hat{k}\}) p_{\hat{k}}(\{\hat{k}\}) d\{\hat{k}\}. \qquad (20) $$

The ensemble average $\mu_\epsilon$ of the error can then be derived using nested conditional averages applied to the two sets of concern in (18) as follows

$$ \mu_{\epsilon} = E[\epsilon] = E_{\mathbf{k}}[E_{w}[\epsilon]]. \qquad (21) $$

From expressing the average at the inner level, it follows that

$$ \begin{aligned} E_w[\epsilon] &= \frac{1}{N} E_w[\mathbf{w}^H \mathbf{\Lambda} \mathbf{w}] - 1 \\ &= \frac{1}{N} \sum_{i=1}^{N} \lambda_i E[|w_i|^2] - 1. \end{aligned} \qquad (22) $$

Reminding that $\mathbf{X}$ is a unitary matrix yields

$$ E[|w_i|^2] = E[|v_i|^2] = 1 \quad \forall i \in [1, N]. \qquad (23) $$

Consequently,

$$ E_w[\epsilon] = \frac{1}{N} \sum_{i=1}^{N} \lambda_i - 1. \qquad (24) $$

The sum over the $\lambda_i$ is equal to $\mathrm{Tr}(\mathbf{P})$, the trace of $\mathbf{P}$. Besides, the main diagonal of $\mathbf{P}$ being made of ones, it follows that

$$ \mathrm{Tr}(\mathbf{P}) = \sum_{i=1}^{N} \lambda_i = N, \qquad (25) $$

which implies in a clear way that

$$ E[\epsilon] = 0. \qquad (26) $$

We can see that the estimator $\hat{\mu}_{p_{\omega}^2}$ is still unbiased for a finite number of plane waves, provided that these latter are sufficiently stirred to justify of a statistical approach that assumes the independence between the sets of random variables. The unbiased nature of the spatial estimator provides a necessary but not sufficient condition to consider spatial-averaging and ensemble-averaging measures as mean-ergodic processes$^{22}$. Indeed, one needs to study the variance of the error to assert the equivalence properties in a more general frame.
### C. Variance of the error

The mean of the error being null according to (26), the variance $\sigma_{\epsilon}^2$ of the error reduces to $E[\epsilon^2]$. Adopting much the same approach as previously, the variance can be expressed using nested conditional averages as

$$ \sigma_{\epsilon}^2 = E_{\mathbf{k}}[E_w[\epsilon^2]]. \quad (27) $$

From (18) we express $\epsilon^2$ as follows

$$ \begin{aligned} \epsilon^2 &= \left( \frac{1}{N} \mathbf{w}^H \mathbf{\Lambda} \mathbf{w} - 1 \right)^2 \\ &= \left( \frac{1}{N} \sum_{i=1}^{N} \lambda_i |w_i|^2 \right)^2 + 1 - \frac{2}{N} \sum_{i=1}^{N} \lambda_i |w_i|^2. \end{aligned} \quad (28) $$

At the inner level of (27) the conditional average on the $\{w\}$ set can be developed as

$$ \begin{aligned} E_w[\epsilon^2] &= E_w \left[ \left( \frac{1}{N} \sum_{i=1}^{N} \lambda_i |w_i|^2 \right)^2 \right] + 1 \\ &\quad - \frac{2}{N} \sum_{i=1}^{N} \lambda_i E_w \left[ |w_i|^2 \right]. \end{aligned} \quad (29) $$

Recalling that the components $w_i$ are normalized random variables, the use of (25) yields the simpler form

$$ E_w[\epsilon^2] = \frac{1}{N^2} E_w \left[ \left( \sum_{i=1}^{N} \lambda_i |w_i|^2 \right)^2 \right] - 1. \quad (30) $$

To obtain a more explicit expression, we expand the term between square brackets as

$$ \left( \sum_{i=1}^{N} \lambda_i |w_i|^2 \right)^2 = \sum_{i=1}^{N} \lambda_i^2 |w_i|^4 + \sum_{i \neq j} \lambda_i \lambda_j |w_i|^2 |w_j|^2. \quad (31) $$

This allows recasting (30) as

$$ E_w[\epsilon^2] = \frac{1}{N^2} \left( \nu_4 \sum_{i=1}^{N} \lambda_i^2 + \nu_2^2 \sum_{i \neq j} \lambda_i \lambda_j \right) - 1, \quad (32) $$

where the $\nu_n$ terms are the following conditional moments defined as,

$$ \nu_n = E_w[|w_i|^n] = E[|w_i|^n]. \quad (33) $$

In order to make (32) more tractable, we use algebraic properties stating that

$$ \sum_{i=1}^{N} \lambda_i^2 = \mathrm{Tr}(\mathbf{P}^T \mathbf{P}) = \|\mathbf{P}\|_F^2 = \sum_{i=1}^{N} \sum_{j=1}^{N} P_{ij}^2, \quad (34) $$

where $\|\mathbf{P}\|_F$ stands for the Frobenius norm of matrix $\mathbf{P}$ Taking the square of (25)

$$ (\mathrm{Tr}(\mathbf{P}))^2 = \sum_{i=1}^{N} \lambda_i^2 + \sum_{i \neq j} \lambda_i \lambda_j = N^2, \quad (35) $$

and using (34), we can show that

$$ \sum_{i \neq j} \lambda_i \lambda_j = N^2 - \|\mathbf{P}\|_F^2. \quad (36) $$

Whence,

$$ E_w[\epsilon^2] = \frac{\|\mathbf{P}\|_F^2}{N^2} (\nu_4 - \nu_2^2). \quad (37) $$

In order to express the variance of the error, we still need to compute the conditional average on the $\{\hat{k}\}$ set, leading to

$$ \begin{aligned} \sigma_{\epsilon}^2 &= E_{\mathbf{k}} \left[ \frac{\|\mathbf{P}\|_F^2}{N^2} (\nu_4 - \nu_2^2) \right] \\ &= \frac{1}{N^2} (\nu_4 - \nu_2^2) E_{\mathbf{k}} \left[ \|\mathbf{P}\|_F^2 \right]. \end{aligned} \quad (38) $$

We can note that the conditional moment on $\{\hat{k}\}$ only applies on the square Frobenius norm since this latter is the only term in (37) depending on the plane-wave directions. Reminding the structure of $\mathbf{P}$, we can restate (34) as

$$ \|\mathbf{P}\|_F^2 = \sum_{i,j} P_{ij}^2 = N + \sum_{i \neq j} P_{ij}^2 \quad (39) $$

The non-diagonal elements $P_{ij}$, denoted by a function $\Phi(\hat{k}_i - \hat{k}_j)$, inherit the iid property from the $\{\hat{k}\}$ random set. Accordingly, we can expand the conditional average $E_{\mathbf{k}}[\|\mathbf{P}\|_F^2]$ as

$$ E_{\mathbf{k}}[\|\mathbf{P}\|_F^2] = N + N(N-1)E_{\mathbf{k}}[\Phi^2(\hat{k}_i - \hat{k}_j)]. \quad (40) $$

Finally, the variance of the error can be expressed as

$$ \sigma_{\epsilon}^2 = \left[ \frac{N-1}{N} G_{\Omega}(\mathbf{a}) + \frac{1}{N} \right] (\nu_4 - \nu_2^2) \quad (41) $$

with,

$$ G_{\Omega}(\mathbf{a}) = \iint_{\Omega} \Phi^2(\hat{k}_i - \hat{k}_j) p_{\hat{k}}(\hat{k}_i) p_{\hat{k}}(\hat{k}_j) d\hat{k}_i d\hat{k}_j, \quad (42) $$

which is the only term depending on the dimensions of $\Omega$, where $\mathbf{a}$ is a vector gathering the sampling-area dimensions $a_v$ used in (12). Using a Monte-Carlo method and assuming a uniform angular distribution, the evolution of $G_{\Omega}$ is plotted in Fig.1, where a line of length $a$ (dashed line), a square of side $a$ (dotted line) and a cube of side $a$ (solid line) have been considered, respectively.

As shown in Appendix A, the term $\nu_4 - \nu_2^2$ in (41) is the relative-ensemble variance $\sigma_{r,p^2}^2$ of the mean-square pressure, defined as,

$$ \sigma_{r,p^2}^2 = E \left[ \frac{\left( p_{\omega}^2(\mathbf{r}) - \mu_{p_{\omega}^2}(\mathbf{r}) \right)^2}{\mu_{p_{\omega}^2}(\mathbf{r})} \right]. \quad (43) $$

If no spatial averaging is performed, then $G_{\Omega}(0) = 1$. It follows that the variance $\sigma_{\epsilon}^2$ of the relative error reduces logically to $\sigma_{r,p^2}^2$, i.e., to the relative variance of
the mean-square pressure in a single point. The relative-ensemble variance has been widely studied in the literature using wave theory$^{7,23}$ and modal theory$^{3-5}$ models and will not be analyzed in detail herein. All these models have been validated experimentally, and show that the relative-ensemble variance $\sigma_{r,p^2}^2$ is inversely proportional to the number of overlapping modes $M_s$ and tends to one for an ideal diffuse field, i.e., for an infinite number of plane waves.

Concerning spatial averaging, the modal models of Davy$^4$ and Weaver$^5$ included the influence of multiple independent receivers, i.e., of independent spatial samples, in the overall relative variance. Accordingly, if we consider a single fixed source and $N_R$ independent receivers, the mean estimated value $\hat{\mu}_{p^2_\omega}$ of $p^2_\omega(r)$ can be deduced using an arithmetic mean, i.e., a discrete version of (7), as

$$ \hat{\mu}_{p^2_\omega} = \frac{1}{N_R} \sum_{i=1}^{N_R} p_{\omega,i}^2, \qquad (44) $$

which is the traditional mean-estimator. Using modal models$^{4,5}$, it follows that

$$ \sigma_e^2 = \frac{1}{N_R} \sigma_{r,p^2}^2, \qquad (45) $$

where $1/N_R$ is referred to as the convergence rate. Within the modal approach, it appears that the more (independent) receivers we take the more accurate the estimated mean value will be. This of course unrealistic because a spatial correlation exists; but as explained in section II modal models cannot include easily the effect of spatial field distribution into the moments of spatial-average estimates.

Using the plane-wave model presented herein, (41) can interestingly be restated as follows

$$ \sigma_e^2 = G(\mathbf{a}, N) \sigma_{r,p^2}^2, \qquad (46) $$

where $G(\mathbf{a}, N)$ is the convergence rate of the spatial estimator expressed as

$$ G(\mathbf{a}, N) = \frac{N-1}{N} G_{\Omega}(\mathbf{a}) + G_{PW}(N), \qquad (47) $$

with $G_{PW}(N) = 1/N$.

Note that on a single point, i.e., for $\mathbf{a} = \mathbf{0}$ the convergence rate $G(\mathbf{0}, N) = 1$ as expected, independently of the number $N$ of plane waves.

It turns out that the convergence rate of the spatial estimator obtained in (47) is conditioned by two terms. The first term, as expected, is related to the influence of the sampling term $G_\Omega(\mathbf{a})$ whose inverse corresponds to an equivalent number of independent samples. The number of plane waves has very little influence on the first term of (47). On the contrary, the second term depends only on the number of plane waves.

The influence of the number of plane waves in the spatial-averaging convergence rate has never been predicted to the best of our knowledge by any study.

<figure><img src="image_15.png" /><figcaption>FIG. 1. Evolution of the spatial term $G_\Omega(\mathbf{a})$ in (42). 1D, 2D and 3D sampling domains have been considered, related to a line stretch of length $a$ (dashed line), to a square (dotted line) and a cube of side $a$ (solid line), respectively.</figcaption></figure>

The study of the accuracy of the spatial-average estimator consists then in analyzing the convergence rate given in (47) which can be simplified as,

$$ G(\mathbf{a}, N) \approx G_{\Omega}(\mathbf{a}) + G_{PW}(N), \qquad (48) $$

since the term $(N-1)/N$ can be reasonably neglected.

Accordingly, to study which of the two terms dominates in the overall estimator variance, the knowledge of $\sigma_{r,p^2}^2$ is not needed since this latter is common. Consequently, we shall focus our analysis in the next section, on a comparison of the two terms of (48).

# IV. ACCURACY OF THE SPATIAL ESTIMATOR

## A. Existence of a lower bound

We are interested in comparing the two terms in (48) in order to see if the contribution of the plane waves may be significant in some cases. This is of relevant importance because the presence of the $1/N$ term in the variance expression acts as a lower bound in the spatial-average estimator accuracy. Indeed, the number of plane waves is related to the room configuration, but is independent of the sampling area $\Omega$. Consequently, if the sampling area is larger, the $G_\Omega$ term will decrease, as seen in Fig. 1, but leaving the $1/N$ term unchanged. In practical terms, raising the sampling area will not diminish continuously the degree of uncertainty deduced from the spatial-average-estimator variance.

To stress the significative risk of this limitation, we study in the next subsection practical room cases for which an assessment of the two contributing terms, i.e., of $G_{PW}$ and $G_\Omega$, is performed.
<figure><img src="image_16.png" /><figcaption>FIG. 2. Reverberation times from Ref. 7. Four configurations have been considered: a large room (solid line), a large damped room (dotted line), a small room lightly damped (dashed line) and a small damped room (dashed-dotted line).</figcaption></figure>

## B. Influence of room parameters

In order to assess quantitatively the role of $G_{PW}$, we consider real room configurations. We use the experimental data presented in Ref. 7 for the case of a small room (40 m³) and a large reverberation room (245 m³) with large fixed diffusers. Both rooms are essentially rectangular. For the small room, lightly damped and damped cases are considered. The corresponding Schroeder's frequencies are 500 Hz and 330 Hz, respectively⁷. For the large room very little damping and large damping cases have been considered. The corresponding Schroeder's frequency are of 310 Hz and 200 Hz, respectively⁷. The measurements of the reverberation time $T_{60}$ were performed on a frequency range [100 Hz : 3.2 kHz] and is shown in Fig. 2.

The configuration of a room has a direct impact on the modal topographies but also on the frequency response and consequently on the number of plane waves composing the sound field. In order to fruitfully use the experimental room parameters we need to establish a way to relate them to the number of plane waves in the following way. At a frequency $f$, we may consider $\beta$ plane waves per resonant mode. In order to assess the number $M$ of modes that are found within the average resonant - 3dB bandwidth $B_M$, we use the simplest form of Weyl's approximation²⁴ of the modal density. The number of plane waves involved in the sound field can then be approximated as

$$N = \beta \cdot M = \beta n(f)B_M \approx \beta \frac{4\pi V f^2}{c^3} B_M, \quad (49)$$

where $n(f)$ is the modal density (in modes per Hertz), $V$ is the volume of the room and $c$ the sound speed. Recalling that the reverberation time $T_{60}$ is related to the average bandwidth of a mode $B_M$, such as $B_M = 6.91/\pi \cdot T_{60}$ yields the following expression of the plane-

<figure><img src="image_17.png" /><figcaption>FIG. 3. Comparison of the lower bound $G_{PW}$ as a function of frequency (lower x-axis) (with ten plane waves per mode) with the 2D-square spatial term $G_\Omega$ as a function of $a/\lambda$ (upper x-axis) for (a) the small room configurations and (b) the large room configurations of Ref. 7. A value of 3.3 m was considered for the parameter $a$.</figcaption></figure>

wave term in (48)

$$G_{PW} = \frac{c^3 T_{60}}{27.64 \beta V f^2}, \quad (50)$$

providing a relation allowing to assess the lower-bound term using the room reverberation time $T_{60}$ and the room volume $V$. For the four room configurations, Fig. 3 shows the evolution of $G_{PW}$ as a function of frequency (lower x-axis), and the evolution of the spatial $G_\Omega$ term as a function of $a/\lambda$ in the 2D case for $3 < a/\lambda < 10$ (upper x-axis) where a value of 3.3 meters has been assumed for the parameter $a$. As already said, the number of plane waves per mode is hard to estimate. The rooms being essentially rectangular we can expect approximately eight plane waves per mode. In order to consider an upper bound we use a value of ten plane waves per mode. We stress that this value being used for all the modes on average, the $G_{PW}$ term may not be regarded as overvalued.

For the case of the small room Fig. 3a shows that the plane-wave term dominates even for frequencies above Schroeder's frequency, whereas for the large room cases described in Fig. 3b, the two terms are of the same order of magnitude. These results can be regarded as a further proof to the well-known fact that the acoustics of small rooms are in a way more complicated than those of large ones. In the present case, it is an illustration that, at a given frequency, the smaller the room the worse the accuracy of the spatial estimate will be; not only because of the relative-ensemble variance as commonly known, but
<figure><img src="image_18.png" /><figcaption>FIG. 4. $G_{PW}$ (bold line) and $G_{\Omega}$ (solid line with circles) as a function of the acoustical volume $V_{\lambda}$ of a room. $G_{\Omega}$ is related to the 2D-section of the sampling volume for which two cases are considered : the sampling volume fits the entire room volume ($\alpha = 1$) and only 10% of the room volume is used in the second case ($\alpha = 0.1$). For $G_{PW}$, frequencies below and above Schroeder's frequencies have been chosen (100Hz and 1kHz, respectively) and three values of reverberation times given in Fig. 2 were used, namely the minimum value of 1 s, the maximum value of about 6 s and an in-between value of 3 s.</figcaption></figure>

also due to the spatial-averaging effect assessed herein. In both cases of Fig. 3, the finite number of plane waves has an influence on the overall convergence rate not only for frequencies near Schroeder's limit, as intuitively expectable, but in fact for a frequency range whose upper limit would be approximately ten times Schroeder's frequency.

At this stage a comparison between $G_{PW}$ and $G_{\Omega}(a)$ has not been performed on the basis of a common parameter. To this end, we consider that the 3D-sampling zone, referred to as $\Omega_{3D}$, is a portion $\alpha$ of the room volume $V$, such as $\Omega_{3D} = \alpha V$. In order to include the influence of the frequency (or equivalently the wavelength $\lambda$) at which the measurement is performed, we introduce the acoustical sampling volume defined as $\Omega_{\lambda}^{3D} = \Omega_{3D}/\lambda^3$. Equivalently, we define the acoustical volume of the room defined as $V_{\lambda} = V/\lambda^3$, yielding

$$G_{PW} = \frac{T_{60} f}{27.64 \beta V_{\lambda}}. \qquad (51)$$

The computation of $G_{\Omega}$ requires the dimensions of the rooms. These latter being not known, we consider the following geometrical ratios 3 : 4 : 5, related to the width $a_1$, the height $a_2$ and the length $a_3$ of the sampling volume, respectively. These ratios account for the absence in practice of degenerate eigenfrequencies.

For the sake of plotting clarity, we restrict the analysis to the $a_2a_3$ section, referred to as $\Omega_{2D}$, of the sampling volume $\Omega_{3D}$. To compute the resulting $G_{\Omega}$ term, two cases of $\Omega_{3D}$ are considered : the case for which the sampling volume fits the total room volume, i.e., for $\alpha = 1$, and the case for which the sampling volume repre-

sents 10% of the room volume, i.e., for $\alpha = 0.1$. Accordingly, Fig. 4 shows the resulting 2D-$G_{\Omega}$ term for $\alpha = 1$ and $\alpha = 0.1$, and the plane-wave term. For this latter, frequencies below and above Schroeder's frequency were considered, 100 Hz and 1 kHz, respectively, and three values of reverberation times were chosen on the basis of the experimental data provided by Fig. 2. The selected values correspond to the minimum of 1 s, to the maximum of about 6 s and to an in-between value of 3 s well suited to the large damped room and small lightly damped room cases.

We can see in Fig. 4 that at a fixed frequency (1kHz in the present case), the decrease of the reverberation time induces a lower convergence rate as stated by (51). If the frequency $f$ varies, so does the acoustical volume since this latter is a function proportional to the cube of $f$. If we refer to the large and small rooms cases, taken at 1 kHz, respective acoustical volumes of 6000 and 991 are found. Using Fig. 4 for $\alpha = 0.1$, we find that $G_{PW} = G_{\Omega} = 10^{-2}$ for the small room; for the large room, $G_{PW} = 2 \ 10^{-3}$ and $G_{\Omega} = 3.5 \ 10^{-3}$. This illustrates the influence of the plane-wave term, since it contributes, for the small room case, to 50% of the overall convergence rate.

<table><caption>TABLE I. For sampling volumes corresponding to 10% ($\alpha = 0.1$) and 100% ($\alpha = 1$) of the small- and large-room volumes, respectively, relative margin of errors $m_{r,\epsilon}$ and the corresponding 2D sampling areas $\Omega_{2D}$ are reported. The number $N_{eq}$ of equivalent independent spatial samples is also reported for $\alpha = 0.1$, assumed as a more realistic case.</caption><thead><tr><th>Configuration</th><th></th><th>$m_{r,\epsilon}$</th><th>$\Omega_{2D}$</th><th>$N_{eq}$</th></tr></thead><tbody><tr><td rowspan="2">small room</td><td>$\alpha = 1$</td><td>20%</td><td>31m<sup>2</sup></td><td>-</td></tr><tr><td>$\alpha = 0.1$</td><td>28%</td><td>6.8m<sup>2</sup></td><td>100</td></tr><tr><td rowspan="2">large room</td><td>$\alpha = 1$</td><td>10%</td><td>105m<sup>2</sup></td><td>-</td></tr><tr><td>$\alpha = 0.1$</td><td>15%</td><td>22m<sup>2</sup></td><td>250</td></tr></tbody></table>

In practice, the spatial-average estimate expressed in (7) can be regarded as a random Gaussian variable whose relative standard deviation $\sigma_{\epsilon}$ conditions the confidence interval, i.e., the accuracy of the estimator. We define the margin of error $m_{\epsilon}$ value as the half-width of the confidence interval. If we consider that 95% of the cases fall into this latter, the margin of error is approximately $m_{\epsilon} = 2\sigma_{\epsilon}$.

Our concern dealing with the comparison of the two terms of (48), the estimator efficiency can be described by the following relative margin of error $m_{r,\epsilon} = m_{\epsilon}/\sigma_{r,p}^2 = 2\sqrt{G_{\Omega} + G_{PW}}$, plotted in Fig. 5 as a function of the acoustical volume of the room, for the 1D (triangles), 2D (circles) and 3D (squares) cases. For the plane-wave term a reverberation time of 3 s at 1 kHz has been used. We consider the cases for $\alpha = 1$ (solid line) and $\alpha = 0.1$ (dashed line) in order to assess the influence of the sampling area on the uncertainty range related to the spatial estimator efficiency.
<figure><img src="image_19.png" /><figcaption>FIG. 5. The relative margin of error related to (48) for a 95% confidence interval and a reverberation time of about 3 s at 1 kHz. Sampling volumes in 1D (triangles), 2D (circles) and 3D (squares) are considered. The case for which the sampling volume $\Omega$ fits the entire room volume ($\alpha=1$) is plotted in solid line - the case for which the sampling volume is 10% of the room volume is plotted in dashed line.</figcaption></figure>

If we refer to the practical cases described previously for which $G_{PW} = G_{\Omega} = 10^{-2}$ for the small room and $G_{PW} = 2 \ 10^{-3}$ and $G_{\Omega} = 3.5 \ 10^{-3}$ for the large room, we find the respective relative margins of error of 20% and 10% for $\alpha = 1$, and 28% and 15% for $\alpha = 0.1$.

Regarding the $G_{\Omega}$ term as an inverse number $N_{eq}$ of equivalent independent discrete spatial samples, 100 and 250 equivalent samples are found for the small and large room, respectively, for $\alpha = 0.1$. All these values are summed up in Tab. I where the equivalent sampling 2D areas $\Omega_{2D}$ have also been reported providing the order of magnitude of the surface related to the relative margins of error previously mentioned.

Interestingly, if the $\alpha = 0.1$ case of the spatial estimator were replaced by an ensemble estimator (using a rotating vane or moving the source and receiver positions as in Ref.7) based on $N_s$ independent samples, we would need, on the basis of the same margins of error, i.e., 28% and 15%, 51 and 177 samples for the small and large room, respectively. For the small room case, to compensate the plane-wave term $G_{PW}$, a spatial-averaging process needs twice more independent samples than an ensemble estimator. It is worth stressing that in the case of an ensemble estimator, each realization corresponds to the case of $G(0, N) = 1$, i.e., the number of plane waves has no influence.

The method presented herein can be then extended in a general way to any configuration. It provides a practical theoretical tool for predictions and estimation procedures applied to the mean-square pressure in room acoustics.

In the general model derived in section III, we did not make any assumption on the degree of diffusion of the field unlike the works dedicated to spatial averaging in the literature$^{13,15-18}$. Indeed, these latter have been carried out within the frame of the ideal diffuse field hypothesis.

To the best of our knowledge, no other studies have been presented in a more general case. Accordingly, it is interesting to analyze how the general expression, given in (46), becomes in the particular case of an ideal diffuse field.

### C. The ideal diffuse field case

#### 1. Relative variance expression

In the case of an ideal diffuse field, i.e., for an infinite number of plane waves, using (47), we can show that the convergence rate $G(\mathbf{a}, N)$ reduces to $G_{\Omega}(\mathbf{a})$. The relative-error variance of the spatial estimator given in (46) can then be expressed as

$$ \lim_{N \to \infty} \sigma_{\epsilon}^{2} = G_{\Omega}(\mathbf{a}) \sigma_{r,p^{2}}^{2}. \qquad (52) $$

The relative-ensemble variance of the mean square pressure $\sigma_{r,p^2}^2$ can be analyzed in a straightforward way. Indeed, the complex plane-waves weights $\{w\}$ in the new basis result from an infinite sum of the $\{\gamma\}$ random variables. Accordingly, using the central limit theorem, we can state that the real and imaginary parts of the set $\{w\}$ follow a normal distribution. Consequently, the elements $|w_i|^2$ follow a chi-square distribution with two degrees of freedom whose moments are such that

$$ \frac{E[|w_i|^4]}{E[|w_i|^2]^2} = \frac{\nu_4}{\nu_2^2} = 2. \qquad (53) $$

Hence, recalling that $E[|w_i|^2] = 1$, the relative variance $\sigma_{r,p^2}^2$ reduces to unity as commonly accepted for the ideal diffuse field case.

We can then express (52) in the following simpler form

$$ \lim_{N \to \infty} \sigma_{\epsilon}^{2} = G_{\Omega}(\mathbf{a}), \qquad (54) $$

which is the expression given by Schroeder$^{13}$. In Schroeder's work$^{13}$ a Gaussian assumption is made on the $\{\gamma\}$ set distribution. Note that, in the model derived in section III such assumption is not necessary to recover the expression of the error variance in the case of an ideal diffuse field.

#### 2. Estimators equivalence

In order to analyze if a spatial- and ensemble-average estimators are equivalent, one needs to compare their convergence rates.

To this end, we introduce the ensemble-mean-estimated value of the mean-square pressure $\bar{\mu}_{p^2}$ defined as,

$$ \bar{\mu}_{p^2} = \frac{1}{M_r} \sum_{i=1}^{M_r} p_{\omega, i}^2, \qquad (55) $$

where $p_{\omega, i}^2$ is the mean-square pressure corresponding to the i-th realization of the ensemble statistics and $M_r$ is
the number of independent realizations. A realization can either be a position of a rotating vane or a pair of source/receiver positions. If the total number $M_r$ of realizations are independent, $\bar{\mu}_{p_\omega^2}$ is a random variable whose centered and normalized form, as done in (13), follows a normal distribution of zero mean and of variance $\sigma_{r,p^2}^2/M_r$.

The spatial-average (centered and normalized) estimate $\epsilon$ follows also a normal distribution with a zero mean. However, the variance and/or convergence rate is not identical in all the spatial-sampling configurations as shown in Fig. 1. Indeed, we observe different slopes preventing a general equivalence principle between spatial- and ensemble-average estimates as the confidence interval is not identical. Accordingly, the claim in Ref. 17 that the spatial-average estimator provides the "true estimate" is not correct in a strict sense although the 1D/2D/3D-spatial-variances $G_\Omega$ tend asymptotically to zero.

In order to carry on a more accurate comparison of the convergence rates of $\hat{\mu}_{p_\omega^2}$ and $\bar{\mu}_{p_\omega^2}$, we focus on the different slopes of $G_\Omega$ for the 1D, 2D and 3D cases which turn to be 0.5, 0.4 and 0.26 per octave, respectively. For the 1D case, the slope is in accordance with the asymptotic expansion of the expression of Ref. 13. For the 2D-square and 3D-cube cases no analytical expressions being available, we have estimated the following asymptotic expressions

$$G_{\Omega}(a) \approx \begin{cases} 0.5/(a/\lambda) & \text{in the 1D case} \\ 0.2133/(a/\lambda)^{1.8} & \text{in the 2D case} \\ 0.0740(a/\lambda)^{-2} & \text{in the 3D case} \end{cases} \quad (56)$$

valid for $a/\lambda > 2$.

These expressions, derived for the first time for the 2D and 3D cases to the best of our knowledge, allow to see the impact of the dimensionality of $\Omega$ on $G_\Omega$ slopes.

The difference of slopes, linked to the difference of exponent power value in the asymptotic expressions of $G_\Omega$, is induced by the spatial correlation of the field in the room. The correlation effects on the spatial variance have already been studied in the cases of a line, a circle and a disk sampling areas$^{15-17}$. For the mean-square pressure, the spatial autocorrelation function $R(k_0\Delta r)$ was shown to be$^{17,25}$,

$$R(k_0\Delta r) = \left( \frac{\sin(k_0\Delta r)}{k_0\Delta r} \right)^2, \quad (57)$$

where $\Delta r$ is the distance between two points.

Using (57), a zero correlation is found for distances of $\lambda/2$, or entire multiples. On the basis of this property, the ensuing idea of using a discrete, rather than a continuous, spatial-average estimator has been studied by Waterhouse$^{18}$ who was interested in studying the different possible arrangements in 1D, 2D and 3D to have uncorrelated discrete samples. It can be shown that the 1D equivalent discrete averaging estimator, i.e., sampling along a line on positions apart from $\lambda/2$ (or multiples), is the only case for which uncorrelated samples can be obtained. Indeed, in the 2D and 3D cases no geometrical configuration can provide a zero correlation between

<figure><img src="image_20.png" /><figcaption>FIG. 6. Variance of the spatial-average estimator in the 1D (triangles), 2D (circles) and 3D (squares) cases with respect to the ensemble-mean estimator (solid line with no markers), as a function of the discrete number of samples.</figcaption></figure>

each pair of points. Consequently, in the 1D case, doubling the number of (independent) samples, i.e. doubling $a/\lambda$, will reduce the estimator variance $\sigma_\epsilon^2$ of a factor 2, confirming the slope obtained in the 1D case in Fig. 1.

To compare $G_\Omega$ with the convergence rate $1/M_r$ of the ensemble estimator, we need to use a discrete spatial-average estimator in order to have the number samples as a common parameter.

It has been shown$^{18}$ that the relative-error variance of the discrete spatial-average estimator can be fairly approximated by $G_\Omega$, as long as the geometric parameter $a$ is greater than $3/2$ $\lambda$. Accordingly, assuming a sufficient number of samples allows using (56) to describe properly the asymptotic convergence rate of the discrete spatial-average estimator. In order to compare the convergence rates of the two estimators, we need to introduce a number of spatial samples whose arrangement will impact inevitably on the accuracy of the spatial-average estimate.

We choose a regular arrangement of the samples separated by a distance $\lambda/2$. This arrangement allows, on the one hand, to have samples easily placed and on the other hand, to have a zero-correlation in many cases. If $N_s$ is the number of samples in the 1D case ensuing from the regular arrangement, $N_s^2$ and $N_s^3$ will be the samples obtained in the 2D-square and 3D-cube cases, respectively. Moreover, if a length $a$ is assumed for the 1D case, the integer part of $1 + a/(\lambda/2)$ correspond to the number of samples. Equivalently, for a sufficient value of $a$, we have approximately $2a/\lambda$ samples along a line or along the side of a square and a cube of length $a$. Consequently, using the asymptotic approximations of $G_\Omega$ we can compute the following asymptotic convergence rates

$$G_{\Omega}(N_s) \propto \begin{cases} N_s^{-1} & \text{in the 1D case} \\ N_s^{-0.9} & \text{in the 2D case} \\ N_s^{-0.666} & \text{in the 3D case} \end{cases} \quad (58)$$
Using the results of $G_\Omega$ shown in Fig.1, we plot in Fig.6 the variance of the relative error of both estimators as a function of the number of samples.

The 1D case and the ensemble estimator being based on independent samples, they converge with the same rate. Whereas, the 2D and 3D cases are affected by spatial correlation. Consequently, their convergence rates are lower. Accordingly, for the scenarios where the ideal diffuse field model can be assumed as reliable, the often-invoked equivalence assumption between spatial-average (in the 2D and 3D cases) and ensemble-average estimators does not hold.

## V. CONCLUSION

The present study showed that the spatial- and ensemble-average estimators cannot be regarded as equivalent, since their convergence rates differ. The results were obtained using a wave-theory approach based on a finite number of plane waves; the plane-wave model is equivalent to a modal expansion, while allowing to explicitly express the spatial dependence of the sound-field distribution, well-suited to the study of the spatial estimator.

For an infinite number of plane waves the spatial-averaging estimates coincide with the results presented forty years ago using an ideal-diffuse field model. In the general case, the present work points out theoretically that the number of plane waves induce a limitation of the spatial-estimator accuracy, emphasizing the uselessness of over-sampling. This effect has never been predicted to the best of our knowledge. On the basis of practical cases, the results show that this effect impacts the convergence rate for frequencies about ten time Schroeder's frequency.

Although the model is physically justified and theoretically reliable, it would be interesting to establish this limitation on the basis of experimental measurements. These studies would provide an experimental evidence of the maximum accuracy that can be expected.

## APPENDIX A: IDENTIFICATION OF THE VARIANCE

In order to identify the relative-ensemble variance $\sigma_{r,p^2}$ of the mean-square pressure defined as,

$$ \sigma_{r,p^2}^2 = E \left[ \left( \frac{p_{\omega}^2(\mathbf{r}) - \mu_{p_{\omega}}^2(\mathbf{r})}{\mu_{p_{\omega}}^2(\mathbf{r})} \right)^2 \right] = E[\epsilon^2], \quad (A1) $$

we must consider the error expression in (13) in a single point. The matrix $P$ is then made of ones. Consequently, the matrix being of unitary rank, we have only a single non-zero eigenvalue $\lambda_1$ such as,

$$ \lambda_1 = N. \quad (A2) $$

Accordingly, (32) can be expressed with the following simpler form,

$$ E_w[\epsilon^2] = \nu_4 - 1. \quad (A3) $$

Recalling that $\nu_2^2 = 1$, the final variance is obtained computing the conditional average on the $\{\hat{k}\}$ set. However, $\nu_4$ is independent of the $\{\hat{k}\}$ set, finally yielding

$$ \sigma_{r,p^2}^2 = \nu_4 - \nu_2^2. \quad (A4) $$

1. M. Schroeder, “Statistical parameters of the frequency response curves of large rooms”, J. Audio Eng. Soc 35, 299–305 (1987).

2. M. Schroeder and K. Kuttruff, “On frequency response curves in rooms. comparison of experimental, theoretical, and monte carlo results for the average frequency spacing between maxima”, J. Acoust. Soc. America 34, 76 (1962).

3. R. Lyon, “Statistical analysis of power injection and response in structures and rooms”, J. Acoust. Soc. America 45, 545 (1969).

4. J. Davy, “The relative variance of the transmission function of a reverberation room”, Journal of Sound and Vibration 77, 455–479 (1981).

5. R. Weaver, “On the ensemble variance of reverberation room transmission functions, the effect of spectral rigidity”, Journal of Sound and Vibration 130, 487–491 (1989).

6. R. Waterhouse, “Statistical properties of reverberant soundfields”, J. Acoust. Soc. America 35, 1894 (1963).

7. F. Jacobsen and A. Molares, “The ensemble variance of pure-tone measurements in reverberation rooms”, J. Acoust. Soc. America 127, 233 (2010).

8. C. Balachandran, “Random sound field in reverberation chambers”, J. Acoust. Soc. America 31, 1319 (1959).

9. F. Jacobsen and T. Roisin, “The coherence of reverberant sound fields”, J. Acoust. Soc. America 108, 204 (2000).

10. H. Kuttruff, *Room acoustics*, chapter 3, 86–88 and 76, 5th edition (Spon Press, Oxon) (2009).

11. W. Joyce, “Sabine’s reverberation time and ergodic auditoriums”, J. Acoust. Soc. America 58, 643 (1975).

12. M. Witting, *Modelling of diffuse sound field excitations and dynamic response analysis of lightweight structures*, chapter 2, 10 (Herbert Utz Verlag, München) (1999).

13. M. Schroeder, “Spatial averaging in a diffuse sound field and the equivalent number of independent measurements”, J. Acoust. Soc. America 46, 534 (1969).

14. M. Schroeder, “Effect of frequency and space averaging on the transmission responses of multimode media”, J. Acoust. Soc. America 46, 277 (1969).

15. D. Lubman, “Spatial averaging in a diffuse sound field”, J. Acoust. Soc. America 46, 532 (1969).

16. D. Lubman, “Spatial averaging in sound-power measurements”, J. Acoust. Soc. America 45, 337 (1969).

17. D. Lubman, R. Waterhouse, and C. Chien, “Effectiveness of continuous spatial averaging in a diffuse sound field”, J. Acoust. Soc. America 53, 650 (1973).

18. R. Waterhouse and D. Lubman, “Discrete versus continuous space averaging in a reverberant sound field”, J. Acoust. Soc. America 48, 1 (1970).

19. P. M. Morse and K. U. Ingard, *Theoretical Acoustics*, chapter 9, 556–558 (McGraw-Hill, Inc., New York) (1968).

20. R. Lyon, “Needed: a new definition of diffusion”, J. Acoust. Soc. America 56, 1300 (1974).

21. P. Chinnery and V. Humphrey, “Experimental visualization of acoustic resonances within a stadium-shaped cavity”, Physical Review E 53, 272 (1996).

22. A. Papoulis, *Probability, random variables, and stochastic processes*, chapter 9, 245–247, 2nd edition (McGraw-Hill, Inc., New York) (1984).

23. F. Jacobsen and A. Molares, “Sound power emitted by a
pure-tone source in a reverberation room", J. Acoust. Soc. America **126**, 676 (2009).

<sup>24</sup> H. Weyl, “Das asymptotische Verteilungsgesetz der Eigenwerte linearer partieller Differentialgleichungen (mit einer Anwendung auf die Theorie der Hohlraumstrahlung)(Asymptotic distribution of the partial differential equations eigenvalues (with an application to black

body theory))”, Mathematische Annalen **71**, 441–479 (1912).

<sup>25</sup> W. Chu, “Note on the independent sampling of mean-square pressure in reverberant sound fields”, J. Acoust. Soc. America **72**, 196–199 (1982).
# Goodness-Of-Fit tests in Reverberation Chambers : Is Sample Independence Necessary?

Florian Monsef*, Ramiro Serra†, Andrea Cozza*

*GeePs UMR 8507, Univ Paris-Sud, CentraleSupelec, CNRS, UPMC, Gif-sur-Yvette, France,
E-mail: florian.monsef@geeps.centralesupelec.fr

†Eindhoven University of Technology, Eindhoven, The Netherlands.

***Abstract***—Goodness-of-fit (GoF) tests are interesting tools to test the overmodedness of the field inside a reverberation chamber (RC). In the definition of GoF tests independent and identically-distributed samples are required. In the present paper we focus on the effect of partial data correlation on GOF tests. We analyze the sensitivity of the most common GOF tests used in electromagnetic compatibility. The results presented herein are based on numerical and experimental data. We show that most of GOF tests are insensitive to correlations up to rather high values depending on the type of GOF test chosen. Implications for practical applications are discussed.

***Index Terms***—Autocorrelation Coefficient, Data Independence, Random Fields, Reverberation Chambers, Statistical Electromagnetism.

## I. INTRODUCTION

Statistical dispersion of field-related quantities in a reverberation chamber (RC) is conditioned by the degree of homogeneity of the field inside the chamber. Ideally, the excitation of a minimum number of modes [1] allows the field to approximatively follow a Gaussian law.

It is common to quantify stirrer efficiency by assessing the number of stirrer steps for which the field can be regarded as reasonably uncorrelated between two contiguous steps [2]. The question of how to estimate this number of positions has been the topic of many papers [3]–[6] and will not be considered here. If the case of the standards are considered, a large number of samples is used in order to obtain a more accurate estimate of the autocorrelation coefficient, from which a decimation factor can then be deduced. This decimation factor allows reducing the number of stirrer positions to the one that ideally corresponds to independent test configurations within the RC, avoiding redundant test results and ultimately saving time.

In order to check the Gaussian fit of the field, goodness-of-fit (GOF) tests can be used as explained in [2]. To the best of our knowledge, it is not clear if the set of data to use must be the decimated one, i.e., the set made up of uncorrelated data, or the complete set of data which provides more samples, even though partially correlated. It is indeed commonly thought that GoF need independent samples to be applied. In practice, samples being correlated, the condition of using independent samples is cut short and one may wonder about the possible use of correlated (and whence dependent) data.

The present work aims to give some insight into this question. To that end, we will recall in Section II the most common

GOF tests used in the context of RCs, while Section III will analyze GOF outcomes according to the autocorrelation model used to generate correlated data. Section III-C will also focus on GOF results but in the case of measured data. Finally, a tentative interpretation of the results will be proposed in Sec. IV.

Our findings show to what extent different GOF tests usually implemented within the EMC community are sensitive to data autocorrelation.

## II. GOODNESS-OF-FIT TESTS

In order to determine if a set of samples belong to a given probability distribution, GOF tests can be used as a statistical decision tool. These tests can be performed, e.g., with the real and/or imaginary parts of the electric field in order to check if these follow a Gaussian law. But these tests can also be applied to field modulus samples in order to check whether these follow a Rayleigh distribution. In both cases, a *null hypothesis*, referred to as $H_0$, is formulated and tested.

There are several existing GOF tests, and, with the same set of data, the outcome can differ from one GOF test to another. This is due to the fact that each GOF is based on a different metric. The question raised herein is not to know which GOF is the more adequate to RC applications but rather to use the most popular ones and to analyze their sensitivity to correlated data. Although not exhaustive, the most commonly encountered ones are the following:

* the $\chi^2$ test [7],

* The Lilliefors (L) test [8],

* the Kolmogorov-Smirnov (KS) test [9],

* the Anderson-Darling (AD) test [5], [10].

The $\chi^2$ test determines whether there is a significant difference between the theoretical frequencies of occurrence and the observed frequencies. It tells if the inevitable difference between the expected and observed values is due to stochastic dispersion in the sample, or if there is a significant underlying difference. The Lilliefors and KS tests both estimate the maximum distance between the empirical cumulative distribution function (cdf) and the expected one. The distance is compared to a critical value, often referred to as $D_{KS}$, in order to decide whether the sample data belong to the expected distribution. The Lilliefors test uses sample moments unlike KS for which the moments of the expected distribution need to be specified *a priori*. The AD test integrates a weighted
difference between the sample cdf and the expected one. The weighted function is such that the test focuses more on the tails of the distribution than in other tests.

In what follows, these GOF tests will be performed on synthetic data generated numerically and on data obtained experimentally.

### III. SENSITIVITY OF GOF ON SYNTHESIZED DATA

In this section we will present an approach where the first lag correlation coefficient is controlled in order to analyze if the outcomes of GOF tests are somewhat related to the degree of correlation.

#### A. Correlation Models

The idea is to generate correlated data to emulate, on average, what one could obtain by using a stirrer in practice. We are interested in generating two different correlation models. This allows us to see if the sensitivity is only linked to the lag-1 correlation, referred to as $ρ_l$, or if the correlation model itself has an impact on GOF outcomes.

We consider a set $x_i$ whose samples follow a standardized normal law. For such a case, let us recall that the lag-1 coefficient is such that $ρ_l = E[x_i x_{i-1}]$, where $E[·]$ is the ensemble average operator.

The first model is a first-order autoregressive (AR(1)) [11] model based on the following recursive equation:

$$ 
\begin{aligned}
x_1 &= \eta_1, \\
x_i &= \rho x_{i-1} + \eta_i \sqrt{1 - \rho^2},
\end{aligned} 
\quad (1) 
$$

where the random variable $η$ is standardized, i.e., $η \sim N(0, 1)$. For such an AR(1) model, we report in Fig. 1 the autocorrelation function obtained (by using 5000 samples) for $ρ = 0.3$ and $ρ = 0.7$, respectively. It is clear that $ρ_l$ corresponds to $ρ$.

The second model considered herein is an AR(2) model for the correlation that may exist between a field value obtained at a given position and those obtained at the two preceding positions. The recursive equation reads in that case,

$$ x_i = \rho x_{i-1} + \rho' x_{i-2} + \gamma \eta_i. \quad (2) $$

Note that (2) can be regarded as a dynamic system that may be unstable<fnref n="1" /> if $\{ρ, ρ'\}$ are not properly chosen. As in (1), the coefficient weighting $η_i$, here $γ$, is set to ensure $x_i$ to be standardized such that $γ = \sqrt{1 - (ρ + ρ')^2}$. This condition also provides a restriction on the possible values of $ρ$ and $ρ'$ since $(ρ + ρ')^2$ must be less than 1. The (somewhat arbitrary) choice made to fulfill the constraints on $ρ$ and $ρ'$ is such that $ρ' = 0.25 ρ$. Although not intuitive, for an AR(2) model, $ρ_l$ is a function of $ρ$ and $ρ'$ such that $ρ_l = ρ/(1 - ρ')$, i.e., is larger than for an AR(1) model for a given $ρ$ and $ρ' \neq 0$. In what follows, the number of data per set will be referred to as $N$.

#### B. Results

By considering a AR(1) model and by using (1), we generated data used as input to the GOF tests. The simulation

<figure><img src="image_21.png" /><figcaption>Fig. 1. Autocorrelation function obtained by using the model given by (1) for $ρ = 0.3$ (asterisks) and $ρ = 0.7$ (squares), respectively.</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>Fig. 2. Rejection rates ($α = 5%$), obtained with an AR(1) model, for different samples sizes and GOF tests as a function of data autocorrelation.</figcaption></figure>

study is carried out by generating sets of $x_i$ with $N = 50$ and $N = 500$ and for 50 values of $ρ$ uniformly distributed over [0:1]. To reduce statistical uncertainties, each realization was repeated a thousand times.

Fig. 2 summarizes the percentage of rejected null hypotheses of normality for different GOF tests under study and different sample numbers, for a significance level $α = 5\%$ (p-value). We stress the fact that these rejection rates are solely coming from data correlation, since the tested data sets come from a normal parent distribution, as described by the process in (1).

It can be seen that the KS test shows a higher sensitivity to data autocorrelation than other GOF tests. It seems also more robust against different sample sizes. From these results, it appears that KS test could be suggested as a way of identifying residual autocorrelation on data, even for low values of $ρ$.

In order to highlight the case of the KS test sensitivity, critical values are plotted in Fig. 3 for a varying number of independent realizations<fnref n="2" /> and different correlation levels. We can see that $D_{KS}$ values are shifted upwards as the correlation rises. For $ρ < 0.5$, we notice that $D_{KS}$ is almost insensitive to the degree of correlation, i.e., is almost similar to the iid case.

Another interesting way to tackle this analysis is to wonder which correlation level corresponds to an average p-value of 5%, the average being performed over sets of p-values taken from tests repeated a thousand times. Fig. 4 shows the results

<footnote n="1">¹if the poles modulus is larger than 1.</footnote>

<footnote n="2">²that correspond to stirrer states ideally independent in practice.</footnote>
<figure><img src="image_23.png" /><figcaption>Fig. 3. The KS threshold obtained for different correlation levels as a function of the number of stirrer states. The threshold obtained for the classic (uncorrelated) case (solid line) is shown as well as the thresholds obtained for values of $\rho$ reported in the figure and running from 0 to 0.9.</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>Fig. 4. Correlation threshold leading to an average empirical p-value of 5% as a function of the number of stirrer states.</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>Fig. 5. Rejection rates obtained with an AR(2) model for $N = 500$. $\chi^2$ (dashed line) AD (dashed dotted line) and Lilliefors (Green) do not seem to be sensitive to correlation, unlike the KS test, for which the AR(1) result has been reported for comparison.</figcaption></figure>

<figure><img src="image_26.png" /><figcaption>Fig. 6. Average lag-1 autocorrelation coefficient of measured data inside an RC at 1GHz, for different levels of sample decimation (different sample size). Vertical bars stand for 95% confidence intervals.</figcaption></figure>

obtained for the different GOF tests (see legend). We observe that the correlation level for non-KS tests is quite high and almost independent of the number of stirrer states.

Now, let us analyze the rejection rate obtained when using an AR(2) model. Fig. 5 shows the corresponding results. It can be seen that non-KS tests (dashed line, dashed-dotted line and green line) are still insensitive to correlation, whereas the KS rejection rate increases faster (solid line) than what was obtained with an AR(1) model (dotted line). In a way, the AR(2) provides a higher average correlation level<fnref n="2" /> which could explain the difference of slopes in the rejection-rate curves obtained for the two AR models considered here.

### C. Experimental Setup and Measurements

We will check now if the sensitivity of the GOF tests studied herein is also observed when experimental data are used.

Measurements were performed in the RC facility of the Eindhoven University of Technology. The chamber is 4 × 5 × 2.85 (h) m with a lowest useable frequency of 200 MHz, according to the criteria in [2]. A signal generator, a RF power

amplifier and a log-periodic antenna (80 MHz - 1 GHz) were used in order to inject electromagnetic energy into the RC at a fixed frequency of 1 GHz. The three cartesian components of the electric field magnitude were measured by an electric field probe located inside the working volume for a total of 3600 stirrer positions (one field sample every 0.1°). As is well known, ideally each field component is expected to follow a Rayleigh distribution. By taking such small stirrer step angles, we intend to generate highly correlated data that are also reasonably Rayleigh distributed as well.

The procedure is to start with the complete set of 3600 samples, estimate $\rho_l$ and apply the different GOF tests. Subsequently, we will decimate the set, in order to get a lower $\rho_l$ (as well as a lower number of samples) and again apply the different GOF tests. We perform this process with increasing decimation until reaching a set of 12 samples. Fig. 6 shows the estimated $\rho_l$ by using the three field components as a function of the sample size. As expected, the degree of correlation lowers as the decimation factor increases.

Three different GOF tests were applied to the increasingly decimated data: the $\chi^2$, the KS, and the AD tests. The Lilliefors test could be applied to test against a Rayleigh

<footnote n="2">²as the for the lag-1 case, for a given $\rho$, the lag-2 coefficient is larger for the case of an AR(2) model than for an AR(1) model; it equals $\rho' + \rho^2/(1 - \rho')$ for the former and only $\rho^2$ for the latter.</footnote>
<figure><img src="image_27.png" /><figcaption>Fig. 7. Rejection rates ($\alpha = 5\%$) for different GOF tests (see legend) as a function of lag-1 autocorrelation for measured data.</figcaption></figure>

distribution but tables needed for that specific case are not available. Fig. 7 highlights the rejection rates as a function of lag-1 correlation for the measured (and decimated) sets of data. Similar conclusions as from the simulation study can be drawn: GOF tests other than the KS appear insensitive to data autocorrelation, even for relatively high values of $\rho_l$, whereas KS tests seem to be as sensitive as what was observed with numerical results.

## IV. DISCUSSION

To understand the pertinence of the results presented in the previous sections, it is worth recalling that the main effect of correlation is to reduce the range value spanned by the successive data samples, by introducing a direct dependence, roughly said, a sort of drag or inertia.

This has a impact on the number of bins that need to be set for the $\chi^2$ test. As long as the correlation allows the data to span a sufficiently wide range of values to balance the frequencies in the different classes according to Cochran's criterion, the $\chi^2$ test remains operational. When the correlation becomes too strong, the values spread collapses to a narrow value interval that makes possible only the definition of a single bin. And then the test fails.

The same idea can be transposed to the AD test where the empirical cdf fails to be representative of the whole real distribution; only a small portion of the cdf becomes accessible and the comparison to a reference distribution becomes no longer possible.

Concerning the KS test, the outcome is conditioned by Kolmogorov statistics. This statistics are sometimes called the Brownian bridge that assumes a random Brownian motion. Variations are due to trajectories, assumed to be independent. This independence is not the same concept required for independent data samples, as repeatedly used throughout the present work. However, a causality seems to exist. The correlation has an impact on the KS critical value, but to the best of our knowledge, there is no simple relation between the proper $D_{KS}$ value and $\rho_l$. From the results presented in the previous section, the typical critical value $1.36/\sqrt{N}$ (solid line in Fig. 3) seems too be low when correlation intervenes, such that $H_0$ is often rejected.

Finally, the Lilliefors test seems to be more robust by using sample moments than the KS test, but it is not easy to be more specific on the reasons that make it more robust.

## V. CONCLUSIONS

Uncorrelated data were shown not to be necessary when performing $\chi^2$, Anderson-Darling, or Lilliefors GOF tests. It seems that these tests fail when a correlation of about 90% is reached. However, the Kolmogorov-Smirnov test was shown to be very sensitive to correlation. These results were obtained by using synthetic data and experimental data. The rejection rate obtained highlighted the same clear trends.

As stated in standards, GOF evaluates the mismatch between the measured and the theoretical (ideal) distribution functions of field or power density. As highlighted in the present work, an RC user can perform some of these GOF tests without caring about the possible correlation that may exist between samples. In other terms, if the test is well chosen, there is no need to decimate the data to have an insight on the overmodedness of the chamber.

As a perspective, the KS GOF test could be thought as a way to detect correlation in scenarios where the data distribution is known in advance.

## REFERENCES

[1] F. Monsef, "Why a Reverberation Chamber Works at Low Modal Overlap," *IEEE Transactions on Electromagnetic Compatibility*, vol. 54, no. 6, pp. 1314-1317, 2012.

[2] *Reverberation chamber test methods*, International Electrotechnical Commission (IEC), Std. 61 000-4-21, 2011.

[3] F. Moglie and V. Primiani, "Analysis of the independent positions of reverberation chamber stirrers as a function of their operating conditions," *IEEE Transactions on Electromagnetic Compatibility*, vol. 53, pp. 288-295, 2011.

[4] S. Pfennig and H. Krauthauser, "A general method for determining the number of independent stirrer positions in reverberation chambers," in *Electromagnetic Compatibility (EMC EUROPE), 2012 International Symposium on*. IEEE, 2012, pp. 1-6.

[5] C. Lemoine, P. Besnier, and M. Drissi, "Estimating the effective sample size to select independent measurements in a reverberation chamber," *IEEE Transactions on Electromagnetic Compatibility*, vol. 50, no. 2, pp. 227-236, 2008.

[6] O. Lundén and M. Backstrom, "Stirrer Efficiency in FOA Reverberation Chambers. Evaluation of Correlation Coefficients and chi-squared tests," in *IEEE International Symposium on Electromagnetic Compatibility*, vol. 1. IEEE, 2000, pp. 11-16.

[7] J. Kostas and B. Boverie, "Statistical Model for a Mode-Stirred Chamber," *IEEE Transactions on Electromagnetic Compatibility*, vol. 33, no. 4, pp. 366-370, 1991.

[8] F. Monsef and A. Cozza, "Goodness-of-fit tests in radiated susceptibility tests," in *Workshop on Aerospace EMC, 2012 Proceedings ESA*. IEEE, 2012, pp. 1-5.

[9] P. Corona, G. Ferrara, and M. Migliaccio, "Reverberating chambers as sources of stochastic electromagnetic fields," *IEEE Transactions on Electromagnetic Compatibility*, vol. 38, no. 3, pp. 348-356, 1996.

[10] V. Mariani Primiani and F. Moglie, "Numerical simulation of reverberation chamber parameters affecting the received power statistics," *IEEE Transactions on Electromagnetic Compatibility*, vol. 54, no. 3, pp. 522-532, 2012.

[11] C. Lemoine, P. Besnier, and M. Drissi, "Using autoregressive models to estimate the number of independent samples available with stirrer rotation in reverberation chamber," in *Electromagnetic Compatibility, 2007. EMC Zurich 2007. 18th International Zurich Symposium on*. IEEE, 2007, pp. 469-472.
# Why a Reverberation Chamber Works at Low Modal Overlap

Florian Monsef, Member, IEEE

**Abstract**—The physical parameter well adapted to assess the degree of overmodedness of a reverberation chamber (RC) is the number $M_M$ of modes overlapping in a mode bandwidth. The lowest usable frequency of an RC often corresponds to a low modal overlap of one or two modes. Notwithstanding, in spite of this poor number of modes the RC still works. We show, using Monte Carlo simulation, that the number of modes must, in fact, not be restrained to $M_M$ and that the number of modes contributing to the field statistics can be, even at low modal overlap, somewhat larger than expected.

**Index Terms**—Cavity resonators, electromagnetic compatibility (EMC), modal analysis, parametric statistics, reverberation chamber (RC).

## I. INTRODUCTION

IN RECENT years, there has been a considerable interest in understanding the physics of reverberation chambers (RC) used for electromagnetic (EM) test methods [1]. Due to the complexity of the system, a statistical approach is often adopted [2]–[5] and considers the EM field as a superposition of normal modes [6] or, alternatively, as a superposition of plane waves [7]. The former allows for a more physical description, whereas the latter can be regarded as a proper first-approximation model only if the RC is well overmoded. For this regime, the commonly accepted statistical distribution of the electric energy follows a chi-squared distribution with 6 degrees of freedom (DOF). This distribution is asymptotic and assumes an infinite number of modes to be uniformly excited at working frequency in a continuous wave mode.

In reality, the frequency response of each mode corresponds to a Lorentzian whose influence becomes negligible beyond a certain bandwidth. Accordingly, as shown in Fig. 1, one may look for the effective bandwidth $B_e$ that includes the dominant modes.

Intuitively, the –3-dB bandwidth of a mode, referred to as $B_M$, can be thought of as being the right frequency bandwidth into which the dominant modes fall. This intuition may probably ensue from the analogy that one may make, consciously or unconsciously, with filtering theory, where the –3-dB bandwidth is regarded as the frequency band beyond which the attenuation is often regarded as sufficient.

Moreover, $B_M$ is an easy assessable quantity in the context of an RC, and the modal overlap, or degree of overmodedness,

Manuscript received September 7, 2012; accepted September 13, 2012. Date of publication October 25, 2012; date of current version December 14, 2012.

The author is with the Laboratoire des Signaux et Systèmes, Département de Recherche en Electromagnétisme, UMR 8506 Univ Paris-Sud, SUPELEC, CNRS 91190 Gif-sur-Yvette, France (e-mail: florian.monsef@lss.supelec.fr).

Color versions of one or more of the figures in this paper are available online at http://ieeexplore.ieee.org.

Digital Object Identifier 10.1109/TEMC.2012.2222028

<figure><img src="image_28.png" /><figcaption>Fig. 1. Illustration of the modal overlap. $B_e$ is the effective frequency extension to assess, centered around a working frequency (vertical solid line).</figcaption></figure>

conveniently refers to the number of modes found within this bandwidth.

The aim of this short paper is to check whether this notion of *modal overlap* suffices to provide the right order of magnitude on the number of effective modes contributing to the field statistics.

The determination of such a number obviously cannot be controlled by experimental means, making us resort to numerical methods with which the number of modes can easily and arbitrarily be set.

In the first part of this paper, the field model adopted will be presented. The second step consists in validating the numerical approach by checking that the statistical moments fit those obtained experimentally. This validation allows to proceed to a parametric analysis aimed at finding the number of excited modes.

## II. FIELD MODEL

To model the physics in an RC, a modal approach is used, as in [4]. The electric field $\mathbf{E}$ at a position $\mathbf{r}$ and at a frequency $f$ can be approximated as follows:

$$ \mathbf{E}(\mathbf{r}, f) \approx \sum_{i=1}^{M_e} \alpha_i \psi_i(f) \mathbf{e}_i(\mathbf{r}), \quad (1) $$

where $M_e$ is the number of modes included in $B_e$, $\mathbf{e}_i$ is the eigenvector of the $i$th mode that is assumed not to have a specific form, $\psi_i$ is the frequency response of the $i$th mode, and $\alpha_i$ is the modal weight, i.e., the coupling constant of the source to the $i$th eigenmode.

Limiting our analysis to a generic position $\mathbf{r}$, the following factorization for the modal topographies can be used

$$ \mathbf{e}_i(\mathbf{r}) = e_i(\mathbf{r})\hat{\xi}_i(\mathbf{r}), \quad (2) $$

where $\hat{\xi}_i(\mathbf{r})$ is a unitary polarization vector assumed to be uniformly distributed over $4\pi$ sr.
Accordingly, (1) can be recast as

$$ \mathbf{E}(\mathbf{r}, f) \approx \sum_{i=1}^{M_e} \tilde{\gamma}_i \psi_i(f) \hat{\xi}_i(\mathbf{r}), \quad (3) $$

where the complex quantities $\tilde{\gamma}_i$ are regarded as equivalent modal weights, whose real and imaginary parts follow a Normal law, respectively. As stated in [4], the law assumption is not constraining; a uniform law would not greatly change the final results.

It is worth recalling that the number of modes intervening in (1) reads

$$ M_e = \int_{f-B_e/2}^{f+B_e/2} m(\phi) d\phi, \quad (4) $$

where $m(f)$ is the modal density in Hz$^{-1}$. A usual approximation of this latter reads

$$ m(f) \approx \frac{8\pi V f^2}{c^3}, \quad (5) $$

where $V$ and $c$ are the volume of the RC and the speed of light, respectively. It is worth stressing that the modal density is, in a strict way, a fluctuating quantity whose mean is given by (5).

### III. NUMERICAL MODEL VALIDATION

The first step before proceeding to a parametric study consists in validating the numerical model based on Monte Carlo (MC) simulations. The effective bandwidth (see Fig. 1) is chosen such as $B_e = k_M B_M$ where $k_M$ is referred to as the relative bandwidth extension and is set to 51. We will see in the next section that this value is large enough, ensuring all dominant modes take part in the field statistics. The validation criterion is based on a comparison between the relative variance $\zeta_W^2$ obtained by numerical simulations and the one obtained experimentally.

#### A. Experimental Setup

The experimental measurements are similar to those used in [4], where both empty and loaded (hybrid absorber made up of four pyramids about 30 cm high, standing in the center of the RC) cases are considered. As a brief recall, the setup takes place in a 13.3-m$^3$ RC equipped with a 100-step mechanical stirrer blade of 50 cm wide; its lowest usable frequency (LUF) being around 550 MHz. The relative variance is studied over the frequency range of 0.7–3 GHz. The measured quality factors, as a function of frequency, allow to compute the modal bandwidth often regarded as the correlation bandwidth in frequency stirring techniques [6].

For the empty and loaded scenarios, the maximum modal bandwidths are 150 and 375 kHz, respectively. On the frequency range of interest, 1000 linearly spaced frequency bins are used, inducing a 2.3-MHz frequency space between each point. This frequency space is much larger than the maximum value of the modal bandwidths previously mentioned, allowing to consider each measured point as uncorrelated between one another. The uncertainty can then be further minimized by applying a moving average over five contiguous points, followed by a decimation whose factor equals 5 accordingly. It is worth recalling that

modal bandwidths and the modal density are two parameters that can considerably fluctuate. The moving average followed by the decimation aim at extracting the mean trend, in accordance with the model used herein.

In practice, the modal bandwidth $B_M$ is related to the composite quality factor $Q$ [8] of the chamber and the working frequency $f$, such as $B_M = f/Q$.

To compute the number of modes $M_M$, we would normally need to use (4) where $B_M$ would replace $B_e$. However, the relative bandwidth $B_M/f_i$ being of the order of 0.07%, the integral can be reduced, and $M_M$ reads

$$ M_M = m(f) B_M. \quad (6) $$

#### B. MC Simulation Setup

In order to lead this study, we need to vary the number of overlapping modes by modifying the modal density. As stated by (5), two parameters can be changed for this purpose: the frequency $f$ (like in the experimental approach) or the volume $V$ of the RC. The former requires us to know how the quality factor varies with frequency. Although this could be done using the quality factor values given in [4], it is more convenient to adopt the latter way, i.e., to vary the volume of the RC at fixed working frequency and fixed quality factor; this also has the main advantage of coping with a single modal bandwidth $B_M$ set arbitrarily herein to 1MHz; the volume $V$ will be chosen to obtain the desired modal overlap $M_M$, and $B_e$ has been set large enough (to 51) to include the dominant modes.

Each MC simulation consisted in generating 20 times 50 000 independent random realizations of the electric field described by (3). The resulting estimated variance is averaged over the 20 values.

The statistical distributions of the different quantities are those assumed in [4] as they have been validated experimentally. Accordingly, the real and imaginary parts of the equivalent modal weights were assumed to follow a Normal law. A uniform distribution was assumed for both the polarization of the modes, over $4\pi$ sr, and the resonance frequencies spread around the working frequency; the uniform distribution of the resonance frequencies is justified by the fact that the maximum frequency extension corresponds to a relative bandwidth of only 1.5%.

#### C. Results

In the experimental case, quality factors being known, the number of overlapping modes in the -3 dB modal bandwidth can easily be computed, using (6). The resulting relative variances as a function of $M_M$ are shown in Fig. 2 for the empty (solid gray line) and loaded (dashed blue line) cases, respectively.

For the frequency range used, the number $M_M$ runs from 1 to 13 and from 2.5–35 for the empty and loaded cases, respectively. This explains that the main part of the curve in the empty scenario is rather on the left part of the plot. The ratio of 2.5 found between the minimum values of $M_M$, for the empty and loaded
<figure><img src="image_29.png" /><figcaption>Fig. 2. Relative variance of the electric energy density obtained by MC simulations (full black circles) with the corresponding (vertical) uncertainty bars, and experimentally for the empty (solid gray line) and loaded (blue dashed line) scenario.</figcaption></figure>

cases, is directly linked to the ratio between the maximum value of $B_M$ obtained for these scenarios.

The results obtained by MC simulations have also been reported in Fig. 2 (full circles) for $M_M = 1, 2, 3, 5, 10, 15, 20, 25, 30$, and 35.

It is worth recalling that when experimental variances are deduced from sets made up of 100 samples, an unavoidable uncertainty has to be taken into account. To estimate the latter, MC simulations were used again; this merely consisted in computing 5000 times an estimated relative variance obtained over sets composed of 100 realizations—the standard deviation, referred to as $\sigma_s$, was then derived from the estimated-values ensemble. Recalling the five-points average process performed on the experimental data, the resulting standard deviation reads $\sigma_s/\sqrt{5}$. The estimated variance being normally distributed, and considering the possible fluctuations on a 95% confidence interval, uncertainty bars of twice the resulting standard deviation are reported in Fig. 2 (vertical bars) accordingly.

We can observe a very good agreement between the results obtained numerically and those obtained experimentally. At large modal overlap, i.e., for $M_M = 30$ and 35, the relative variance $\zeta_W^2$ obtained by MC simulation tends to the asymptotic value of $1/3$ predicted by the central-limit theorem.

## IV. EFFECTIVE BANDWIDTH AND NUMBER OF CONTRIBUTING MODES

The aim of this section is to assess more finely the effective bandwidth $B_e$ that has to be considered, if one wants to determine the proper number of excited modes. It is worth recalling that although the number of modes has to be regarded over an infinite frequency range, they do not all contribute in the same manner to the field statistics.

The analysis is based on a parametric study which consists in varying the relative frequency extension $k_M$, or equivalently, the effective bandwidth $B_e$, and this, for different values of $M_M$.

In order to present the different cases on a single plot, we have normalized the relative variance $\zeta_W^2$ to its convergence value, i.e., to the value reported in Fig. 2; this quantity, referred to as

<figure><img src="image_30.png" /><figcaption>Fig. 3. Normalized relative variance for different values of $M_M$ as a function of $k_M$, the frequency extension. The values of $M_M$ reported are 1 (solid line), 2 (blue dashed line), 3 (red dotted line), 5 (dashed-dotted line), 10 (magenta solid line with circles), 15 (green dashed line with squares), 20 (brown solid line with asterisk), 25 (blue dashed line with triangles), 30 (red solid line with diamonds), and 35 (dashed line with crosses). The horizontal bold dashed line stands for an accepted error of 2% on the relative variance $\zeta_W^2$.</figcaption></figure>

$\sigma_N^2$ is shown in Fig. 3 for $M_M = 1$ (solid line), 2 (blue dashed line), 3 (red dotted line), 5 (dashed-dotted line), 10 (magenta solid line with circles), 15 (green dashed line with squares), 20 (brown solid line with asterisks), 25 (blue dashed line with triangles), 30 (red solid line with diamonds), and 35 (dashed line with crosses); all these are plotted as a function of $k_M$ varying from 1 to 51.

For $M_M = 1$, the difference between the values of $\sigma_N^2$ obtained at $k_M = 41$ and $k_M = 51$ corresponds to a shift of 0.086%. This shift can clearly be considered as negligible, and the value obtained at $k_M = 51$ can therefore be regarded as a proper reference. It is noteworthy that the case of $M_M = 1$ has been chosen because it corresponds, as observed in Fig. 3, to the worst case, i.e., to the case that has the slowest convergence rate.

The effective bandwidth $B_e$ can now be determined, but is related, in any case, to an arbitrarily chosen criterion. In this study, we determine the effective bandwidth such that the number $M_e$ of modes included in it, contribute at a 98% level to the field uniformity. In other terms, $B_e$ is chosen with an accepted error of 2% on the relative variance $\zeta_W^2$. This error level of 2% is reported (dashed bold horizontal line) in Fig. 3.

The intersections obtained between each curve and the 2% threshold line provide the minimum value of $k_M$, i.e., the number of modal bandwidths $B_M$ that have to be considered.

Fig. 4 shows the resulting necessary frequency extension to consider as a function of the number $M_M$. The plot shows, quite clearly, that at low modal overlap, the modal bandwidth is not sufficient to properly assess the number of modes taking part to the field statistics. Indeed, for a modal overlap of 1, we can see that more than $16B_M$ have to be considered, confirming that $B_M$ must not be regarded as a filter and must not be restricted to its -3-dB bandwidth.

The resulting number of contributing modes can easily be computed by combining (4) and (6), such as

$$M_e = k_M M_M \quad (7)$$
<figure><img src="image_31.png" /><figcaption>Fig. 4. Minimum value of $k_M$ to consider for an accepted error of 2% on the relative variance $\zeta_W^2$. This value corresponds to the number of modal bandwidth $B_M$ that includes the modes contributing to 98% of the relative variance $\zeta_W^2$.</figcaption></figure>

<figure><img src="image_32.png" /><figcaption>Fig. 5. Number of modes included in the effective bandwidth $B_e$ using (7) (left y-axis). Error made on the relative variance $\zeta_W^2$ if the modal bandwidth $B_M$ is regarded as the right frequency range to use to assess the number of contributing modes.</figcaption></figure>

and the values are reported in Fig. 5 (left y-axis).

We can see that more than 16 modes contribute to 98% of the field uniformity level when a single mode overlaps the modal bandwidth $B_M$.

Although the modes excitation is obviously not uniform, the order of magnitude found for $M_e$ highlights the fact that an acceptable field uniformity can be obtained in an RC.

In the case of our chamber, $M_M = 1$ is obtained at 550 and 700 MHz for the loaded and empty cases, respectively, i.e., for frequencies approaching the LUF. In practice, this shows that the insertion of some losses into the RC can be quite beneficial to

enhance the field uniformity for wanted and unwanted emissions tests.

As shown in Fig. 4, at a higher modal overlap $M_M = 35$, the frequency extension tends to 2, i.e., the modal bandwidth $B_M$ is still not the right bandwidth to consider to assess the number of contributing modes; according to Fig. 5, 35 modes would be considered instead of about 80. This corresponds to more than a factor 2. This factor must not, however, be misleading. Indeed, the quantity of interest is the field uniformity, i.e., the relative variance $\zeta_W^2$.

Accordingly, the pertinent error to look at is the error made on $\zeta_W^2$ whose values are given in Fig. 3 for $k_M = 1$, and summarized in Fig. 5.

For $M_M = 1$, we observe an error of 86% that decreases when the modal overlap rises. We can see that considering $B_M$ at high modal overlap can be more justified; an error of (only) 6% is found for $M_M = 35$.

## V. CONCLUSION

The numerical method presented herein offers a very flexible way to determine the number of modes contributing to the field statistics. It has been shown that for a single mode overlap, the number of contributing modes is more than 15 times larger, highlighting the reasons why an RC can be used at low modal overlap.

This study was extended to larger modal overlap cases, providing a quantitative assessment of the number of modes on a wider range.

## REFERENCES

[1] Reverberation chamber test methods, Int. Electrotech. Commission (IEC), Standard 61 000-4-21, 2003.

[2] L. Arnaut, "Statistics of the quality factor of a rectangular reverberation chamber," *IEEE Trans. Electromagn. Compat.*, vol. 45, no. 1, pp. 61–76, Feb. 2003.

[3] M. Hoijer, "Maximum power available to stress onto the critical component in the equipment under test when performing a radiated susceptibility test in the reverberation chamber," *IEEE Trans. Electromagn. Compat.*, vol. 48, no. 2, pp. 372–384, May 2006.

[4] A. Cozza, "The role of losses in the definition of the overmoded condition for reverberation chambers and their statistics," *IEEE Trans. Electromagn. Compat.*, vol. 53, no. 2, pp. 1–12, May 2011.

[5] C. Lemoine, P. Besnier, and M. Drissi, "Investigation of reverberation chamber measurements through high-power goodness-of-fit tests," *IEEE Trans. Electromagn. Compat.*, vol. 49, no. 4, pp. 745–755, Nov. 2007.

[6] T. Lehman, "A statistical theory of electromagnetic fields in complex cavities," *Interaction Notes*, vol. 494, 1993.

[7] D. Hill, "Plane-wave integral representation for fields in reverberation chambers," *IEEE Trans. Electromagn. Compat.*, vol. 40, no. 3, pp. 209–217, Aug. 1998.

[8] B. Liu, D. Chang, M. Ma, and U. S. N. B. of Standards, *Eigenmodes and the composite quality factor of a reverberating chamber*. Nat. Bureau of Standards, 1983.
Relative variance of the mean-squared pressure in multimode media:
former approaches rehabilitated

Florian Monsef, Andrea Cozza<sup>a)</sup>

Département de Recherche en Électromagnétisme, L2S, UMR8506

Univ Paris-Sud, SUPELEC, CNRS

3 rue Joliot Curie - 91190 Gif-sur-Yvette, France

Dominique Rodrigues,Patrick Cellard,Jean-Noel Durocher

Laboratoire National de métrologie et d'Essais

29, avenue Roger Hennequin - 78197 Trappes cedex

(Dated: August 19, 2014)

The commonly accepted model of the relative variance of transmission functions in room acoustics is the one derived by Weaver which aimed at including the effects of correlation between eigenfrequencies. An analytical expression of the relative variance was advantageously derived by means of an approximated correlation function.The relevance of the approximation used for modeling such correlation is questioned here. Weaver's model was motivated by the fact that earlier models derived by Davy and Lyon, based on independent eigenfrequencies assumptions, led to an overestimation with respect to relative variances found in practice. It is shown that this overestimation is due to an improper definition of the ensemble average applying to the eigenfrequencies that turns out to be at the origin of an inconsistency clearly proven in the present work. Another definition is used that settles the inconsistency and the relative variance values are found to be in good agreement with experimental data. These results rehabilitate former approaches that were based on independence assumptions between eigenfrequencies.

PACS numbers: 43.55.Cs, 43.55.Br, 43.55.Gx

# I. INTRODUCTION

Ensemble variances of the mean-squared pressure and of the sound power have been a topic of interest in room acoustics for more than half a century¹⁻⁶. After the first approach formulated by Schroeder¹,⁷ valid under the diffuse-field regime assumption¹, Lyon⁸ widened the ability of describing the variability of sound-field related quantities for frequencies below Schroeder's frequency.

To that end, he made use of modal theory, which allowed him to introduce the effect of limited modal overlap in a more formal way. His work is, in our opinion, a work of reference where statistical analysis of input impedances and transmission functions have been studied for 2D and 3D structures. His analysis considered the case for which eigenfrequencies are independent, i.e., for which the spacing of adjacent eigenfrequencies follow a Poisson distribution; he also studied the case of a point source on a plate for which the spacing statistics follow a distribution, referred to as the nearest-neighbor (NN) distribution², for which correlation between eigenfrequencies exist. This distribution was improperly regarded some time ago⁹ as a correct model of the level-repulsion effect; for chaotic cavities, the Wigner distribution turns out to be the proper theoretical distribution to consider as stated by random-matrix theory(RMT) [see section 3.1 of Ref.10]. Moreover, it is accepted today that chaotic cavities are systems belonging to the Gaussian

Orthogonal Ensemble (GOE).

Relative-variance models were then extended by Davy² to cover the case when transmission functions are averaged over multiple (independent) source and receiver positions, and this, for both distributions. Concerning the NN distribution, Davy explained that Lyon's formula was not correct and proposed a corrected one².

Some years later, Weaver⁴ proposed another expression of the relative variance of transmission functions. His work was motivated, on the one hand, by the fact that Davy's formula led to an overestimation of relative variances observed in practice, and on the other hand, because the modal-correlation function for GOE systems is rather different from the one obtained by using the NN distribution assumed by Davy² and Lyon⁸. The good agreement found with experimental makes it the accepted model today in room acoustics community. The arguments given by Weaver are all the more justified that the Wigner distribution related to GOE systems was in good agreement with experimental results for mesoscopic systems such as microwave¹¹ and acoustics chaotic cavities¹²,¹³.

In fact, if one analyzes closely the correlation model used by Weaver⁴, it appears that a somewhat crude approximation of the GOE two-point correlation function is used. To highlight this point, we will recall in section II the different correlation models considered in the literature, and this, in order to provide a better understanding of the limits of Weaver's approximation.

Besides, it is worth recalling that Lyon showed (for a 2D case) [pp.560-561 of Ref.8] that the effect of level repulsion is not that relevant at low modal overlap when

<footnote>a)florian.monsef@lss.supelec.fr</footnote>
the variances of input impedances are studied for instance. Because of the high interest for RMT during the two last descales, this finding has been somewhat forgotten apparently. Although based on the NN distribution, the effect of level repulsion highlighted by Lyon motivated the present work for the case of transmission functions. It is worth recalling that the effect of level repulsion shown by Lyon at high modal overlap was challenged by Davy [pp. 458 of Ref.2] but not the one at low modal overlap. This is why Lyon's result raised some doubts and questioning about the overestimation found when comparing the analytical variability expression of transmission functions found at low modal overlap with experimental results, when Poisson distribution is considered. These doubts are also supported by the fact that the relative variance of transmission functions derived under the Poisson assumption has recently led to an analytical expression that was found to be in good agreement for the case of microwave cavities¹⁴,¹⁵.

To highlight these points, the present paper aims at revisiting in section III the derivation of the relative variance of the mean-square pressure by assuming the eigenfrequencies as independent. As we will see, the result differs from Davy's and Lyon's expression. We will present in section IV the origin of the discrepancy and we will show that the revised expression is found to be in good agreement with the experimental results presented in section V.

## II. CORRELATION MODELS

We briefly recall in the present section the concepts needed to understand why, from a theoretical point of view, the claim made by Weaver in the first part of his article⁴ is correct. First, we need to recall what is meant by correlation model. To that end, we refer to Lyon's work⁸ who gave an interesting contribution on that issue and invite the reader to be aware of the aforementioned work.

We adapt Lyon's approach for the case of a "bilateral" eigenfrequencies spectrum regarded as a random process, referred to as $x$, and characterized by a correlation function, referred to as $\Phi_x(\omega)$. When the process $x$ is composed of Dirac's train where each Dirac position is a random variable related to an eigenfrequency, $\Phi_x(\omega)$ can be expressed as,

$$ \Phi_x(\omega) = n_s(\delta(\omega) + S(|\omega|)), \quad (1) $$

where $n_s$ is the modal density. Constant for a 2D cavity, the modal density can be regarded as locally constant for room cases when a single-tone excitation is considered, an assumption also used by Weaver⁴; $S(\omega)$ is the probability density of finding two distinct eigenfrequencies spaced by $\omega$.

If eigenfrequencies are regarded as uncorrelated the function $S(\omega) \equiv S_P(\omega)$ is such that [sec. V of Ref.8],

$$ S_P(\omega) = n_s, \quad (2) $$

and the result for the Poisson case is obtained. The correlation in (2) states that the probability density is uniform

for any spacing - it means that we have the same probability density of finding two eigenfrequencies at the same frequency than at two different frequencies.

For the case of the NN distribution, the function $S(\omega) \equiv S_{NN}(\omega)$ is such that [see Sec. V of Ref.8],

$$ S_{NN}(\omega) = n_s(1 - \exp(-4n_s\omega)). \quad (3) $$

As explained by Weaver, for GOE systems the function $S(\omega)$ to consider is the two-point correlation function, commonly referred to as $R_2(\omega)$, and expressed as,

$$ S_{GOE}(\omega) = R_2(\omega) = n_s(1 - Y_2(\omega)), \quad (4) $$

where $Y_2(\omega)$ is the cluster function whose expression is given by [pp.109 of [10],[4]<fnref n="1" />]

$$ Y_2(\omega) = \left( \frac{\sin(\pi n_s \omega)}{\pi n_s \omega} \right)^2 + \left[ \frac{\pi}{2} \operatorname{sgn}(n_s \omega) - \operatorname{Si}(\pi n_s \omega) \right] \times \left[ \frac{\cos(\pi n_s \omega)}{\pi n_s \omega} - \frac{\sin(\pi n_s \omega)}{(\pi n_s \omega)^2} \right]. \quad (5) $$

It appears clearly that the correlation function is indeed different from the NN case. The important feature of this function is the repulsion effect found at low distance (see pp. 109 of Ref.10) such that the first-order approximation of the two-point correlation function $R_2(\omega)$ reads,

$$ S_{GOE}(\omega) = R_2(\omega) \approx \frac{\pi n_s \omega}{4}, \quad \text{for } n_s \omega \ll 1, \quad (6) $$

that can be compared, by using (3), to the first-order approximation that one would have with the NN distribution,

$$ S_{NN}(\omega) \approx 4n_s\omega, \quad \text{for } n_s\omega \ll 1. \quad (7) $$

The slope in (6) is larger than for the GOE case, which means that the repulsion effect is weaker for the NN case.

Now, if we consider the approximation used by Weaver, the function $S(\omega) \equiv S_W(\omega)$ is such that,

$$ S_W(\omega) \approx n_s(1 - \delta(n_s\omega)). \quad (8) $$

This approximation was presented as valid at high overlap, an argument that we do not consider appropriate for the following reasons. First, if we relate what is meant by high overlap to a diffuse-field regime, the order of magnitude of $n_s$ varies from 3 to 10 according to Schroeder's criterion¹⁶. Second, the approximation could make sense for high distance values $\omega$, but it is in contradiction with the distance range over which the level-repulsion effect is to be expected, i.e., for low values of $\omega$ as stated by Eq. 6 and Eq. 7.

Apart from these arguments, if one compares (8) to (2), both expressions are quite similar except that (8)

<footnote n="1">¹ Note however that the $\pi$ term preceding the sgn(·) operator is missing in [4]</footnote>
has a Dirac's distribution that does not allow one to find two eigenfrequencies at the same value. In other terms, expression (8) is the correlation function that one would obtain with a Poisson-like distribution that does not allow eigenfrequencies degeneracy, i.e., does not allow two eigenfrequencies to be equal.

It appears that the repulsion effect of GOE systems, as the one described by (6), cannot be invoked by using such approximated model.

This is probably the reason why Lyon chose the NN distribution that allowed him to derive an explicit analytical expression of $\Phi_x$. Although the repulsion effect assumed by Lyon is not quantitatively correct, his study on the repulsion effect is of interest and is summarized in Fig. 14 of his work that shows that the variance at low modal overlap of input impedances is not affected by the level repulsion.

Given that Weaver's model is close to a Poisson model and given the study of Lyon on the effect of the repulsion level, the overestimation found for the case of the relative variance of transmission functions in rooms found by assuming independent eigenfrequencies is quite surprising. This motivated to redo the derivation of such quantities in the next section.

### III. RELATIVE VARIANCE OF TRANSMISSION FUNCTIONS

#### A. Sound field

In order to derive the relative-variance expression with the ability to make some connections with Davy's (and Lyon's) works easier, we will start our derivation by adopting the notations used in Ref.2.

For a single-tone excitation, the sound pressure at working angular frequency $\omega_w$, is expressed<fnref n="2" /> as follows,

$$p(\mathbf{x}, \omega_w) = \sum_m -i\omega\rho c^2 \frac{b_m p_m(\mathbf{x})}{\omega_w^2 - \omega_m^2 - 2i\omega_w \gamma_m}, \quad (9)$$

where $\rho$ is the ambient density, $c$ the sound velocity, $p_m(\mathbf{x})$ is the $m^{th}$ eigenmode function, $\gamma_m$ is its related damping constant and $b_m$ is the coupling coefficient of the volume velocity $Q(\mathbf{r})$ on the $m^{th}$ mode at the source location, expressed as,

$$b_m = \frac{1}{V} \int_{(V)} Q(\mathbf{r}) p_m(\mathbf{r}) d\mathbf{r}, \quad (10)$$

where $V$ stands for the effective volume of the source that contributes to the excitation of the $m^{th}$ mode.

Using the usual hypothesis stating that, on the one hand, damping is much smaller than $\omega_w$ and, on the other hand, that the main contribution comes from $\omega_m$ close to angular working frequency $\omega_w$, the sound pressure $p(\mathbf{x}, \omega_w)$ can be approximated as follows [as in Eq. (A23),

Ref.2²],

$$p(\mathbf{x}, \omega_w) \approx \frac{\rho c^2}{2V} \sum_m p_m(\mathbf{x}) p_m(\mathbf{y}) \times (h(\omega_w - \omega_m) - ig(\omega_w - \omega_m)), \quad (11)$$

where,

$$\begin{aligned} h(\omega) &= \frac{\gamma}{\omega^2 + \gamma^2} \\ g(\omega) &= \frac{\omega}{\omega^2 + \gamma^2}. \end{aligned} \quad (12)$$

The rest of the derivation will differ from Davy's and Lyon's approach in the sense that the modal expansion will still be based on an infinite number of modes as required by modal theory.

For brevity, we group the frequency response in a single function, referred to as $\psi(f)$, linked to $h(\omega)$ and $g(\omega)$ such that,

$$\psi(\omega) = h(\omega) - ig(\omega). \quad (13)$$

Moreover we introduce equivalent modal weights defined as,

$$c_m(\mathbf{x}) = p_m(\mathbf{x})b_m, \quad (14)$$

and regarded as centered random variables.

So, at working angular frequency $\omega_w$ the sound pressure can then be restated in a more concise form as follows,

$$p(\mathbf{x}, \omega_w) \approx A \sum_{m=1}^{\infty} c_m(\mathbf{x}) \psi(\omega_w - \omega_m), \quad (15)$$

where $A = \rho c^2/2$.

It is worth recalling that Davy worked separately with the functions $U$ and $W$ defined as,

$$\begin{aligned} U(\omega) &= \mathfrak{Re}(p(\mathbf{x}, \omega_w)/A) \\ W(\omega) &= \mathfrak{Im}(p(\mathbf{x}, \omega_w)/A), \end{aligned} \quad (16)$$

where $\mathfrak{Re}$ and $\mathfrak{Im}$ refers to the imaginary and real part of a complex function, respectively.

#### B. Relative variance derivation

Our aim is to derive the relative variance of the mean-squared pressure defined as,

$$\mathrm{RelVar} (|p(\mathbf{x}, \omega_w)|^2) = \frac{\langle |p(\mathbf{x}, \omega_w)|^4 \rangle}{\langle |p(\mathbf{x}, \omega_w)|^2 \rangle^2} - 1. \quad (17)$$

<footnote n="2">² except that we directly deal here with the sound field and not the transmission function; this one required in former works an assumption of a point-like source which is not assumed herein.</footnote>
In order to have more concise expressions we adopt the following notations for the rest of the derivation,

$$ 
\begin{aligned}
E_{|\psi|^2} &= \left\langle \sum_m |\psi(\omega_w - \omega_m)|^2 \right\rangle \\
E_{|\psi|^4} &= \left\langle \sum_m |\psi(\omega_w - \omega_m)|^4 \right\rangle \\
E_{|\psi|^2|\psi|^2} &= \left\langle \sum_{k \neq l} |\psi(\omega_w - \omega_k)|^2 |\psi(\omega_w - \omega_l)|^2 \right\rangle,
\end{aligned} 
\quad (18) $$

and $\mu_n = \langle |c_m|^n \rangle$, where the brackets $\langle \rangle$ stand for the average over an ensemble of rooms sharing the same characteristics; in practice the ensemble can be taken in a single room into which the source(s) and receiver(s) are displaced so as to modify the boundary conditions and sample spatially the sound field.

We need then to establish the ensemble average of $<|p(\mathbf{x}, \omega)|^2>^2$ and $\langle|p(\mathbf{x}, \omega)|^4\rangle$, respectively. For the former it is straightforward to show that,

$$ \langle |p(\mathbf{x}, \omega_w)|^2 \rangle = \mu_2 E_{|\psi|^2}. \quad (19) $$

For the latter, $|p(\mathbf{x}, \omega_w)|^4$ can be expressed as,

$$ 
\begin{aligned}
|p(\mathbf{x}, \omega_w)|^4 &= \sum_m |c_m|^4 |\psi(\omega_w - \omega_m)|^4 \\
&\quad + 2 \sum_{k \neq l} |c_k|^2 |c_l|^2 |\psi(\omega_w - \omega_k)|^2 |\psi(\omega_w - \omega_l)|^2 \\
&\quad + 2 \sum_m |c_m|^2 |\psi(\omega_w - \omega_m)|^2 \\
&\quad \times \sum_{k \neq l} c_k c_l \psi(\omega_w - \omega_k) \psi^*(\omega_w - \omega_l) + A_{ijkl},
\end{aligned} 
\quad (20) $$

where $A_{ijkl}$ is a sum over four indices with, at least, one of them different from the others.

Given the assumption of independence between the modal parameters and recalling that $\{c_m\}$ sets are made up of centered random variables it follows that,

$$ \langle |p(\mathbf{x}, \omega_w)|^4 \rangle = \mu_4 E_{|\psi|^4} + 2\mu_2^2 E_{|\psi|^2|\psi|^2} \quad (21) $$

The relative variance can be recast at this stage as follows,

$$ \mathrm{RelVar} (|p(\mathbf{x}, \omega_w)|^2) = \frac{\mu_4 E_{|\psi|^4}}{\mu_2^2 E_{|\psi|^2}} + 2 \frac{E_{|\psi|^2|\psi|^2}}{E_{|\psi|^2}} - 1. \quad (22) $$

The critical part deals with the computation of ensemble averages involving the power of $\psi$. To highlight this point let us start with the term $E_{|\psi|^2} = \sum_m |\psi(\omega_w - \omega_m)|^2$ in (19) that corresponds, for a given realization, to the sum of the black dots shown in the upper plot in Fig. 1. Under the assumption of quasi identical $\gamma$ of the most contributing modes, this sum is also obtained approximately by rearranging these dots on the frequency response of the closest mode to working

<figure><img src="image_33.png" /></figure>

<figure><img src="image_34.png" /><figcaption>FIG. 1. The effect of modal overlap at a given angular frequency centered about the working angular frequency $\omega_w$ (upper plot) transposed at the scale of the frequency response of a single mode (lower plot).</figcaption></figure>

frequency, as shown by the lower plot in Fig. 1, so that,

$$ 
\begin{aligned}
E_{|\psi|^2} &\approx \left\langle \int |\psi(\omega - \omega_w)|^2 \sum_m \delta(\omega - \omega_m) \mathrm{d}\omega \right\rangle \\
&= \iint p_{\omega_m}(\omega_m) |\psi(\omega - \omega_w)|^2 \delta(\omega - \omega_m) \mathrm{d}\omega \mathrm{d}\omega_m,
\end{aligned} 
\quad (23) $$

where $p_{\omega_m}(\omega_m)$ is the probability density function (pdf) of the $\omega_m$.

The first line of (23) shows how the density of eigenstates $\sum_m \delta(\omega - \omega_m)$ gets involved in $E_{|\psi|^2}$. The second line allows computing $E_{|\psi|^2}$ provided that the pdf of the $\omega_m$ is known.

To clarify this point let us recall that the ensemble of realizations on which the brackets $\langle \rangle$ apply, are the ensemble of pairs of positions source/receiver in a room. The displacement of the receiver (microphone) and above all of the source (given its larger size) in the room modifies the boundary conditions.

Modifications of these boundary conditions have a direct impact on modal topographies, but also on the eigenfrequencies. On an ensemble of realizations, it is reasonable to consider that the eigenfrequency of a given mode will vary on a bounded interval $\Delta\omega(\omega)$ that can be fre-
quency dependent and that should, at least on average, be the eigenfrequencies mean spacing, i.e., the inverse of the modal density referred to as n whose first-order approximation for 3D structures is such that [pp.78, Ref.17],

$$n(\omega) \approx \frac{V\omega^2}{2\pi^2 c^3}, \qquad (24)$$

in modes/(rad/s).

If we consider that an eigenfrequency $\omega_m$ is uniformly distributed over $\Delta\omega$, it follows that $p_{\omega_m}(\omega_m) = n(\omega_m)$.

Noting that the double integral in (23) can be simplified thanks to the delta function, $E_{|\psi|^2}$ can be recast as,

$$\begin{aligned} E_{|\psi|^2} &\approx \int \frac{1}{\Delta\omega(\omega)} |\psi(\omega - \omega_w)|^2 d\omega \\ &= \int n(\omega) |\psi(\omega - \omega_w)|^2 d\omega \end{aligned} \qquad (25)$$

As already pointed out in Ref.8, this kind of integral does not have a convergence problem since the modal density for the rooms varies as the square of the frequency and $|\psi|^2$ varies asymptotically as $1/f^4$. Given the sharpness of the resonances adopted within our low damping assumption with respect to the angular frequencies of interest, (25) can be further simplified by considering the modal density “sampled” at working frequency, such that,

$$E_{|\psi|^2} \approx n(\omega_w) \int |\psi(\omega - \omega_w)|^2 d\omega = \frac{\pi n(\omega_w)}{\gamma}, \qquad (26)$$

since $\int |\psi(\omega - \omega_w)|^2 d\omega$ corresponds to $h\#h(0) + g\#g(0) = \pi/\gamma$ where $\#$ corresponds to the convolution operator in Davy's notations (see (A37) in Ref.2).

The same reasoning can be followed with $E_{|\psi|^4}$ that reads,

$$E_{|\psi|^4} \approx n(\omega_w) \int |\psi(\omega - \omega_w)|^4 d\omega = \frac{\pi n(\omega_w)}{2\gamma^3}, \qquad (27)$$

since $\int |\psi(\omega - \omega_w)|^4 d\omega$ corresponds to equation (A37), i.e., $(1/\gamma) \times \pi/2\gamma$.

It is easy to show that the ratio $E_{|\psi|^4}/E_{|\psi|^2}$ intervening in (22) is such that,

$$\frac{E_{|\psi|^4}}{E_{|\psi|^2}^2} = \frac{1}{2\pi n(\omega_w)\gamma} = \frac{1}{\pi M} = \frac{1}{M_s}, \qquad (28)$$

where $M_s$ refers to the number of modes overlapping in the statistical bandwidth of a mode, i.e., is $\pi$ times the number M of modes overlapping in the -3 dB modal bandwidth.

Adopting the same approach we expand in appendix the derivation of $E_{|\psi|^2|\psi|^2}$, that turns out to be,

$$E_{|\psi|^2|\psi|^2} \approx E_{|\psi|^4} \left( \frac{E_{|\psi|^2}}{E_{|\psi|^4}} - 1 \right) \qquad (29)$$

Setting $K_c = \mu_4/\mu_2^2$, the relative variance for a single-tone excitation reads,

$$\mathrm{RelVar} (|p(\boldsymbol{x}, \omega_w)|^2) = 1 + \frac{K_c - 2}{M_s}. \qquad (30)$$

The usual K factor being defined as,

$$K = \frac{\langle p_m^4(\boldsymbol{x}) \rangle}{\langle p_m^2(\boldsymbol{x}) \rangle^2} = \frac{\langle p_m^4(\boldsymbol{y}) \rangle}{\langle p_m^2(\boldsymbol{y}) \rangle^2}, \qquad (31)$$

we can relate $K_c$ to the usual K factor such that $K_c = K^2$ and the relative variance reads,

$$\mathrm{RelVar} (|p(\boldsymbol{x}, \omega_w)|^2) = 1 + \frac{K^2 - 2}{M_s}. \qquad (32)$$

Comparing (32) to Davy's (and Lyon's) formula that was found to be$^2$,

$$\mathrm{RelVar} (|p|^2) = 1 + \frac{K^2}{M_s}, \qquad (33)$$

a clear difference appears in the second term of the right-hand term.

This difference is quite surprising since both models have been derived on an independence assumption of the eigenfrequencies. We explain the origin of this discrepancy in the next section.

## IV. ORIGIN OF THE DISCREPANCY

Although both approaches assume a Poisson distribution, unlike in Davy's paper$^2$, the derivation presented herein does not truncate the modal expansion defined in (9) at any stage. This has the major advantage to avoid the issue of defining an arbitrary criterion that will serve to truncate the modal expansion.

This notwithstanding, we need to clarify the origin of the discrepancy. Given the time that has elapsed since the publication of Ref.2 and Ref.8, the reader may be skeptical and surprised that this discrepancy had not been reported since then. In order to raise any possible doubts, an easy and irrefutable way of checking this inconsistency is to perform Monte Carlo (MC) simulations.

These consisted in computing (9). Specific statistical distributions for each modal parameter were assumed and eigenfrequencies were regarded as independent. An ensemble of values of the mean-square pressure was obtained. The corresponding relative variance was then computed and compared to the analytical expressions given by (32) and (33). We assumed $p_m(\boldsymbol{x})$ and $p_m(\boldsymbol{y})$ to be normally distributed; the resulting K factor is therefore equal to 3, which corresponds to the commonly accepted value$^{18}$.

Concerning the MC computations, the question that arises deals with the frequency range over which the eigenfrequencies must be computed. Let us consider a frequency range of about $k_M$ modal bandwidths about a given working frequency. Fig. 2 shows the evolution of