# PARTICIPATION FRANCAISE A LA CAMPAGNE EUROPEENNE D'ETUDE DE LA STRATOSPHERE POLAIRE EASOE

Jean-Pierre POMMEREAU

Service d'aéronomie du CNRS, BP3, 91371 Verrières le Buisson

L'objectif de la campagne EASOE (European Arctic Stratospheric Ozone Experiment) qui s'est déroulée en hiver 1991-92 en Arctique, était l'étude de la photochimie de la stratosphère et en particulier la quantification de la destruction de l'ozone qui peut y prendre place. Conduite avec le soutien de l'Union Européenne et des agences nationales des pays européens, elle a mis en oeuvre une grande variété de systèmes d'observations - 18 stations au sol, 41 vols en ballon, 3 avions de recherche - associés à un large effort de modélisation et d'interprétation. Dans ce contexte, les équipes françaises soutenues par le programme PAMOY réunissant le CNRS, le CNES et les Ministères de l'Environnement et de l'Éducation Nationale, ont largement contribué aussi bien aux efforts d'observation qu'à l'interprétation de résultats. Ceux-ci maintenant publiés dans une édition spéciale du Geophysical Research Letters, ont montré que le potentiel de destruction de l'ozone par les chlores en Arctique était important, mais que la réduction globale de l'ozone en 1992 avait été faible du fait de la date relativement avancée du réchauffement stratosphérique. Cependant, avec l'augmentation prévisible de la teneur en chlore de la stratosphère au cours des dix ans à venir, les résultats de EASOE démontrent que cette réduction pourrait devenir très importante une année où la température de la stratosphère demeurait froide jusqu'à une date avancée dans la saison.

The aim of the European Arctic Stratospheric Experiment (EASOE) which took place in the winter 1991-92 in the Arctic, was the study of the photochemical processes involved in the destruction of stratospheric ozone and the quantification of its global reduction. During the winter. Supported by the European Union and the national agencies of the european countries, the campaign combined a large variety of observational systems - 18 ground stations, 41 balloon flights and 3 research aircraft - and a significant effort for the interpretation and modelisation of the results. In this context, the French groups supported by the national cooperative programme PAMOY of the CNRS, CNES and the Ministères of Environment and Education, contributed significantly to the observations as well as the interpretation of the results. The last now published together in a special issue of the Geophysical Research Letters, have demonstrated the high potential of chlorine for destroying ozone in the Arctic, although the global reduction in the winter 1991-92, was found relatively small because of the early date of the warming of the stratosphere. However, because of the predicted increase of chlorine during the next decade, the results of EASOE demonstrate that the depletion might become very large in a year where the stratosphere would re main cold late in the season.

## 1. INTRODUCTION

Bien qu'aucun trou d'ozone n'ait été observé en Arctique et que par conséquent le risque d'une dégradation rapide de l'ozone stratosphérique de grande envergure soit faible dans l'hémisphère nord, une réduction sensible de l'ozone total aux latitudes moyennes est rapportée à la fois par les observations satellitaires et le réseau de stations au sol de l'Organisation Mondiale de la Météorologie. Elle est maximale en hiver et au printemps avec une amplitude de 0.6 à 0.8 % par an. L'objectif premier de la campagne Européenne EASOE (European Arctic Stratospheric Ozone Experiment) soutenue par la Direction des Technologies de l'Environnement de la Commission Européenne et les agences nationales - CNRS, INSU, CNES, Ministères de l'Environnement et de l'Éducation pour la France - était de tenter de mieux comprendre à la fois le mécanisme et l'amplitude exacte de la réduction d'ozone en Arctique. La campagne d'observations s'est déroulée en hiver 1991-92.

## 2. LES MOYENS

Le projet EASOE a réuni environ 150 chercheurs expérimentalistes et modélistes, principalement Européens mais également Américains, Russes, Japonais et Néo-Zélandais, mobilisés pour de nombreuses observations de terrain pratiquées entre le mois de Novembre 1991 et le mois de Mars 1992 et leur interprétation ultérieure.

Le dispositif d'observation comportait: a) des mesures depuis le sol (lidar, spectromètres uv-visibles, interféromètres infra-rouges, sondeurs micro-ondes et sondages ozone) en dix huit

stations réparties du Groenland à la Sibérie et du Spitzberg à la Grèce; b) la participation de trois avions de recherche, un Transall Allemand emportant un lidar, un spectromètre uv-visible et un interféromètre IR, un Falcon 50 équipé d'un radiomètre micro-ondes pour la mesure de CIO et le Fokker Français ARAT emportant le lidar Léandre pour l'étude des nuages stratosphériques; c) 41 vols en ballons stratosphériques opérés par les équipes du CNES depuis le centre Suédois de l'ESRANGE à Kiruna avec la collaboration du Central Aerological Observatory de Moscou pour la récupération des instruments dans la péninsule de Kola; d) des mesures d'ozone total déduites du sondeur TOVS à bord des satellites météorologiques opérationnels à orbite polaire de la NOAA et produites dans un délai d'une demi-journée par le Centre National de Recherches Météorologiques de Météo-France. Les expériences en ballon consistaient soit en des vols répétitifs d'instruments légers (spectromètre uv-visible pour la mesure de O₂, NO₂ et OCIO, bouteilles de prélèvement de gaz de longue durée sources utilisés comme traceurs du vortex polaire, hygromètre à point de givre) dans des conditions météorologiques variées, soit encore en des vols de nacelles lourdes et longues à préparer, mais permettant la mesure d'espèces difficiles à observer: NOx, HNO₃, ClONO₂, HCl ..

La campagne comme sa préparation étaient conduites par un petit groupe d'organisation appelé "Core Group" en collaboration avec le centre de coordination de Cambridge.

Un élément original de EASOE et qui a largement aidé à l'interprétation ultérieure des observations, a été la mise en place d'une banque de données centrale au Norwegian Institute for Air Research près d'Oslo, dans laquelle les données étaient introduites par chaque groupe et accessibles à tous sous la seule
<table><caption>Tableau 1. Epériences et modèles français associés à EASOE</caption><thead><tr><th>Expérience</th><th>Moyen</th><th>Objectif</th><th>Laboratoire</th><th>Equipe Scientifique</th></tr></thead><tbody><tr><td>SAOZ</td><td>Spectro.uv-vis<br>(4 stations)</td><td>O<sub>3</sub>, NO<sub>2</sub><br>PSC</td><td>SA / CNRS<br>DMI, CAO, NILU</td><td>F. Goutail,<br>J. P. Pommereau,<br>A. Sarkissian</td></tr><tr><td>OHP</td><td>Sol. Lidar</td><td>O<sub>3</sub>, aerosols</td><td>SA / CNRS</td><td>S. Godin, C. David</td></tr><tr><td>ELSA</td><td>Sol. Lidar<br>Sodankyla<br>Finlande</td><td>O<sub>3</sub>, aérosols</td><td>SA / CNRS<br>IROE / Italie<br>U. Berlin</td><td>S. Godin, C. David,<br>C. Flesia, G. Visconti,<br>L. Stefanutti, J. P. Wolf</td></tr><tr><td>SAOZ-B</td><td>Ballon<br>Spectro. uv-vis</td><td>O<sub>3</sub>, NO<sub>2</sub>,<br>OCIO</td><td>SA / CNRS</td><td>J. P. Pommereau<br>J. Piquard</td></tr><tr><td>ELHYSA</td><td>Ballon<br>Hygromètre</td><td>Vapeur d'eau</td><td>LMD</td><td>J. Ovarlez<br>H. Ovarlez</td></tr><tr><td>LOA</td><td>Ballon, Radiomètre<br>Polarimètre</td><td>Aérosols et<br>PSC</td><td>LOA Lille</td><td>M. Herman</td></tr><tr><td>LPMA</td><td>Ballon<br>FTIR</td><td>HCl, HNO<sub>3</sub></td><td>LPMA / Paris 6</td><td>C. Camy-Peyret</td></tr><tr><td> OzONOXY</td><td>Ballon<br>Nacelle<br>Multi-Instr.</td><td>ClONO<sub>2</sub><br>O<sub>3</sub>, NOy<br>HNO<sub>3</sub><br>aérosols<br>CH<sub>4</sub>, CFC, N<sub>2</sub>O</td><td>U. Denver<br>SA / CNRS<br>U. Nagoya<br>U. Denver<br>KFA Julich<br>U. Houston</td><td>F.J. Murcray<br>P. Aimedieu<br>Y. Kondo<br>F.J. Murcray<br>W.Schmidt<br>W.D. Sheldon</td></tr><tr><td>LEANDRE</td><td>Avion ARAT</td><td>PSC<br>Lidar</td><td>SA / CNRS<br>Obs. Neuchatel</td><td>S. Godin, C. David<br>G. Mégie, D. Haner<br>C. Flesia</td></tr><tr><td>TOVS</td><td>NOAA-12</td><td>Ozone total</td><td>CNRM, Météo-France</td><td>D. Cariolle, F. Karcher</td></tr><tr><td> Modèle<br>Lagrangien</td><td>Modélisation<br>Photochimie</td><td>NOx</td><td>SA / CNRS</td><td>E. Lateltin<br>J. P. Pommereau, H. Le<br>Texier</td></tr><tr><td> MCG</td><td>Modélisation</td><td>Dynamique</td><td>LPCE / CNRS<br>ONERA</td><td>M. Pirre, R. Ramaroson</td></tr><tr><td> ARPEGE</td><td>Modélisation</td><td>PSC, O<sub>3</sub></td><td>CNRM<br>Météo-France</td><td>D. Cariolle<br>Y. Orsolini, M. Déqué</td></tr><tr><td> Modèle<br>Chimie-Transport</td><td>Modélisation</td><td>Chlore</td><td>CNRM<br>Météo-France</td><td>D. Cariolle<br>M. Chipperfield P. Simon</td></tr><tr><td> Théorie</td><td>Modélisation</td><td>Ondes de gravité</td><td>LMD</td><td>H. Teitelbaum, J. Ovarlez</td></tr></tbody></table>

réserve de la signature d'un protocole visant simplement à protéger chaque auteur de l'utilisation intempestive de ses propres données.

### 3. LA CONTRIBUTION FRANÇAISE

Les équipes françaises ont largement contribué au succès de EASOE aussi bien sur le plan scientifique que sur le plan opérationnel. Six laboratoires se sont engagés dans des programmes de mesure au sol, en avion ou en ballon ou à partir de satellites, ou encore dans la l'interprétation des données à l'aide de modèles de simulation. La Division Ballons du CNES a assuré l'ensemble des opérations en vol des ballons en Scandinavie, la Division technique de l'INSU, celles du Fokker ARAT basé à l'aéroport de Kiruna et trois représentants français ont participé à la conception, à la mise en place et à la conduite de la campagne au sein du "Core Group".

#### 3.1. Observations au sol.

3.1.1. Le réseau de spectromètres SAOZ. Utilisé pour la première fois en 1988 en Arctique à Sodankyla (Finlande), dans le cadre de la campagne Franco-Allemande d'étude de la stratosphère arctique CHEOPS III, le spectromètre SAOZ a été conçu au laboratoire pour la mesure des colonnes totales d'ozone et de NO₂ tout au long de l'année et par tous les temps à haute latitude. Cette propriété due à l'emploi d'une méthode d'observation de la lumière solaire diffusée au zénith au cours du crépuscule pratiquée pour la première fois dans le domaine visible, a été démontrée auparavant [Pommereau et Goutail, 1988a et 1988b], comme a été démontrée sa capacité à détecter la présence de nuages stratosphériques qui se forment parfois en hiver et encore une fois par tous les temps [Sarkissian et al., 1991].

Quatre stations SAOZ ont été exploitées durant ESAOE: à Sodankyla en Finlande avec le Finnish Meteorological Institute (FMI); à Scoresbysund sur la côte Est du Groenland avec le
Danish Meteorological Institute (DMI); à Zhigansk (123°E, 67°N) en Sibérie Orientale, région on l'on observe en moyenne un maximum d'ozone ainsi qu'un maximum de température stratosphérique à l'opposé de ce qui est rapporté en Arctique Européen, en collaboration avec le Central Aerological Observatory (CAO) de Moscou; et à bord du navire météorologique Norvégien "Polarfront" en station au milieu de l'Atlantique Nord (02° E, 66°N), en amont du dispositif EASOE d'observation de la stratosphère en Scandinavie et en collaboration avec le Norwegian Institute for Air Research (NILU).

Les études de ces données menées par la suite au laboratoire ont porté sur deux sujets: les variations de NO₂ en Arctique et la détection des nuages stratosphériques. Les observations de NO₂ ont démontré une réduction de 30% de la densité de colonne du constituant par rapport aux années antérieures, résultat de la conversion hétérogène des NOx en acide nitrique; l'absence de dénitrification globale du vortex en Arctique durant EASOE comme cela est observé chaque année en Antarctique et enfin une forte et inattendue corrélation entre la colonne de NO₂ et la température de la stratosphère qui suggère soit un mécanisme de saturation de la conversion hétérogène des NOx dépendant de la température, soit un mécanisme de retour à la phase gazeuse non identifié pour l'instant [*Goutail et al., 1994*].

L'analyse d'un indice de couleur du ciel au crépuscule a montré que les premiers signes précurseurs de l'arrivée des aérosols volcaniques du Pinatubo sont apparus en Arctique dès le mois de Septembre 1991 et qu'ils s'étaient peu à peu épaissi en automne aussi bien dans qu'à l'extérieur du vortex polaire dans la basse stratosphère, masquant de possibles PSC en dessous de 22 km. Au dessus de cette altitude, où les aérosols étaient absents, la formation de PSC aurait été observée si elle s'était produite. En réalité, un seul cas de nuage stratosphérique élevé a été observé et uniquement à Sodankyla à la fin du mois de Décembre, contrairement aux plus fréquents événements des années antérieures [*Sarkissian et al., 1994*].

### 3.1.2. Lidars à l'Observatoire de Haute Provence et à Sodankyla.
Le système lidar permettant de mesurer l'ozone stratosphérique de l'OHP ainsi qu'un second, installé au sol à Sodankylä (Finlande, 67°N), dans le cadre du projet coopératif ELSA impliquant 5 groupes européens (Service d'Aéronomie, l'IROE de Florence, l'Observatoire de Neuchâtel, la Freie Universität de Berlin et l'Université de l'Aquila, Italie), ont participé à EASOE [*Godin et al., 1994a*]. Ce projet avait pour objectif l'étude corrélée de la distribution d'ozone stratosphérique et des nuages stratosphériques polaires. Il impliquait la mise en oeuvre sur le même site de 2 systèmes lidar multi-longueurs d'onde permettant de mesurer simultanément l'ozone et les aérosols stratosphériques. Le système installé par le CNRS et émettant à 289, 299, 308 et 355 nm avait pour objectif de mesurer par diffusion élastique différentielle, la distribution verticale de l'ozone de 5 à 40 km. La présence du nuage volcanique du Pinatubo a cependant considérablement perturbé les mesures dans le domaine d'altitude du nuage, situé dans la basse stratosphère précisément là où l'on cherche à mesurer l'effet sur la concentration d'ozone de l'activation des composés chlorés induite par la présence de nuages stratosphériques polaires. L'incertitude sur les mesures lidar correspondant à l'altitude du nuage volcanique (environ 15 - 22 km) est telle qu'il est difficile d'utiliser ces mesures pour estimer une diminution de la concentration d'ozone induite par les nuages stratosphériques polaires ou les aérosols volcaniques dans la basse stratosphère.

L'analyse des observations à la longueur d'onde de référence (355nm) a permis toutefois d'étudier l'évolution des aérosols au-dessus de Sodankylä pendant la campagne EASOE. Le nuage détecté dès le début juillet à l'Observatoire de Haute-Provence, s'est rapidement propagé vers le nord puisqu'il a été observé

début Août à Andoya (69°N, 16°E) et au Spitzberg (79°N). Il était systématiquement présent pendant toute la durée de la campagne EASOE. Une comparaison des mesures d'aérosols stratosphériques effectuées pendant la campagne à différentes latitudes a permis de constater une différence significative dans le domaine d'altitude occupé par les aérosols volcaniques : contrairement aux stations situées en latitude moyenne qui ont pu observer des aérosols jusque vers 700 K (27 km) d'altitude, les stations situées à l'intérieur du vortex polaire n'ont détecté le nuage volcanique que jusque vers 450 K (19 km). La variabilité du rapport de diffusion est liée à la présence ou non du vortex polaire au-dessus de Sodankylä, comme l'indiquent les valeurs de vorticité potentielle superposées sur la figure. Ainsi au début du mois de Janvier, alors que le vortex est localisé au-dessus de la Scandinavie les valeurs mesurées sont plus faibles et l'altitude du nuage plus basse. Ces observations illustrent d'une part l'efficacité des processus de transport entre les régions de latitude moyenne et les régions polaires dans le domaine d'altitude situé entre la tropopause et 450 K et mettent en évidence d'autre part les propriétés isolantes du vortex polaire au-dessus de cette altitude [*Neuber et al., 1994*].

## 3.2. Observations en ballons.

41 vols en ballon ont été effectués durant la campagne, répartis entre deux types d'expériences qui correspondent à deux approches différentes d'observation : le vol répétitif de petits instruments de 20 à 50 kg sous des ballons de faible volume (10 000 à 35 000 m³) dans des conditions météorologiques variées pour tenter de comprendre les causes que régissent les modifications de concentration des espèces et le vol de charges plus complexes et plus lourdes, jusqu'à 400 kg, sous des ballons de 35 000 à 100 000 m³ de volume pour la mesure de la distribution verticale d'espèces difficiles à observer.

Deux instruments français, les spectromètres SAOZ et l'hygromètre Elhysa relèvent de la première approche et ont effectué de nombreux vols, respectivement 11 et 8, alors que trois autres plus complexes ont pu voler 2 à trois fois (Tableau 2).

### 3.2.1. L'instrument SAOZ du Service d'Aéronomie du CNRS est un spectromètre uv-visible destiné à la mesure des profils verticaux de l'ozone, NO₂, O₄ et OCIO et de l'extinction par les aérosols. Il a pu voler 11 fois entre le 28 Novembre et le 19 Mars, à l'intérieur et à l'extérieur du vortex, dans des conditions variées de température. Les résultats de ces vols interprétés à l'aide d'un modèle photochimique Lagrangien, ont permis de quantifier la réduction de la concentration des oxydes d'azote par conversion hétérogène en acide nitrique à la surface des gouttelettes d'aérosols volcaniques [*Pommereau et Piquard, 1994a, Lateltin et al., 1994*]. Ils ont également permis d'obtenir la première observation de la distribution verticale du dioxyde de chlore OCIO résultant de "l'activation" des chloros et des bromes localisées dans la tranche d'altitude la plus froide où la température est suffisamment basse pour permettre la formation de nuages stratosphériques [*Pommereau et Piquard, 1994b*].

### 3.2.2. L'instrument ELHYSA du Laboratoire de Météorologie Dynamique, est un hygromètre à point de givre qui permet des mesures in-situ de la concentration de la vapeur d'eau, paramètre clé dans le calcul de la formation des nuages stratosphériques ou PSC. Les 7 vols effectués entre Novembre et Mars ont démontré une remarquable stabilité de la vapeur d'eau en dehors du vortex, au dessus de 16 km d'altitude (4 à 5 ppmV), mais une large variabilité en dessous, en relation avec l'origine des masses d'air. L'advection d'air tropical s'identifie par exemple par la présence d'air sec. Un autre aspect des résultats obtenus a été l'observation de couches atteignant de très basses températures induites par la propagation d'une onde de relief et indicatives de la formation de
<table><caption>Tableau 2. Expériences françaises en ballon</caption><thead><tr><th>Instrument</th><th>Objectif</th><th>Masse</th><th>Ballon</th><th>Nb de vols</th></tr></thead><tbody><tr><td>SAOZ</td><td>Profils de O<sub>3</sub>, NO<sub>2</sub>, OCIO, O<sub>4</sub></td><td>19 kg</td><td>10 000 m<sup>3</sup></td><td>11</td></tr><tr><td>ELHYSA</td><td>Profils de vapeur d'eau</td><td>50 kg</td><td>35 000 m<sup>3</sup></td><td>7</td></tr><tr><td>LOA</td><td>Aérosols et PSC</td><td>120 kg</td><td>35 000 m<sup>3</sup></td><td>2</td></tr><tr><td>LPMA</td><td>HNO<sub>3</sub>, HCl, ClONO<sub>2</sub></td><td></td><td>100 000 m<sup>3</sup></td><td>2</td></tr><tr><td>OZONOXYPOLY</td><td>O<sub>3</sub>, NO, NOy, HNO<sub>3</sub>, N<sub>2</sub>O, CH<sub>4</sub>, CFC, Aérosols</td><td>450 kg</td><td>100 000 m<sup>3</sup></td><td>3</td></tr></tbody></table>

nuages stratosphériques. D'autre part, il a été montré que la descente diabatique de l'air à l'intérieur du même vortex conduisait à la réhydratation de la stratosphère [Ovarlez et Ovarlez, 1994].

### 3.2.3. L'instrument du Laboratoire d'Optique Atmosphérique de l' Université de Lille est constitué d'un radiomètre et d'un polarimètre à deux longueurs d'onde, qui permet de restituer la distribution verticale de la dimension et de la nature (indice de réfraction) des particules diffusantes, observées au cours de l'ascension du ballon. Il a pu effectuer deux vols, dont l'un en combinaison avec un compteur de particules de l'université du Wyoming, monté sur la même nacelle. Le résultat principal est l'observation à tous niveaux jusqu'à 20 km environ de particules de l'ordre de 0.5 µm de diamètre dont l'indice de réfraction est compatible avec celui d'un mélange à 60% d'eau et d'acide sulfurique. Les caractéristiques des aérosols déduites de la polarisation étaient pratiquement les mêmes pour les deux vols. En terme de rayon effectif pour la diffusion, le meilleur accord correspond à des particules de rayon $r_{eff} = 0.525$ µm et d'indice de réfraction $m=1.37$. L'épaisseur optique verticale des aérosols déduite des mesures est particulièrement élevée; de 0.05 à 0.1 à 1 µm avec une dépendance spectrale de -1.5. De plus, lors d'un vol à basse température, des couches très structurées formées de plus grosses particules indiquant la formation d'un nuage stratosphérique ont été observées vers 19-20 km d'altitude [Herman et al., 1992].

### 3.2.4. L'instrument du Laboratoire de Physique Moléculaire et Applications (LPMA) est un spectromètre infrarouge a transformée de Fourier fonctionnant en absorption et utilisant le soleil comme source. En optimisant l'heure de lancement, on parvient à réaliser des observations à la montée ainsi qu'au plafond, en profitant alors de la géométrie habituelle de l'occultation solaire au coucher. L'instrument est embarqué sur une nacelle pointée approximativement dans la direction du soleil. Elle emporte également un petit spectromètre infrarouge refroidi, CAESR, construit par l'Université de Denver, qui mesure la distribution verticale de HNO<sub>3</sub> à la montée en émission, dans la direction opposée au soleil.

Deux vols ont été effectués. Au cours du premier (5 fev 1992), l'instrument LPMA opérait dans la région 7-12 µm et les les spectres observés confirment les valeurs obtenues par CAESR pour HNO<sub>3</sub> [Murcray et al., 1994].

Au cours du second vol (14 mars 1992), un profil vertical de HCl jusqu'à 26 km a été obtenu en couvrant la région 2.9-3.6 µm.

Ce profil permet, avec celui de ClONO<sub>2</sub> obtenu lors du vol très rapproche de l'instrument MIPAS, de faire un bilan de la partition des chloros dans le vortex juste à la fin de la période froide. On a également mesuré les profils de CH<sub>4</sub> et N<sub>2</sub>O lors du deuxième vol qui a vu également la perte de la nacelle à la suite de la rupture d'un des éléments de la chaîne de vol et d'une chute libre de 30 km.

### 3.2.5. La nacelle multi-instrumentée OZONOXYPOLY du Service d'Aéronomie est conçue comme une plate-forme pouvant emporter différents instruments suivant les objectifs poursuivis: un ozonomètre du Service d'Aéronomie, un instrument de mesure de NO et NOy par chimiluminescence de l'Université de Nagoya au Japon, un compteur d'aérosols du même laboratoire, le spectromètre infra-rouge de l'Université de Denver déjà cité, et des bouteilles de prélèvement du KFA Julich pour la mesure des espèces sources. Trois vols ont été effectués dans des combinaisons instrumentales différentes pour l'observation de NO et des aérosols [Aimedieu et al., 1994]; du rapport entre les oxydes d'azote NOy et leur source N<sub>2</sub>O, sensible au transport diabatique vertical et à une dénitrification éventuelle - absente dans ce cas - [Kondo et al., 1994]; ou encore de la distribution de HNO<sub>3</sub> [Murcray et al., 1994].

## 3.3. Expérience Avion

Le lidar Léandre du Service d'Aéronomie embarqué sur le Fokker 27 français ARAT avait pour objectif l'observation de la distribution des couches de diffusion stratosphériques à moyenne échelle. L'utilisation de deux canaux de polarisation croisée à la longueur d'onde 532 nm, permet une discrimination entre les particules quasi-sphériques des aérosols volcaniques et les particules glacées non-sphériques des nuages stratosphériques. L'appareil a effectué deux campagnes en Scandinavie du 2 au 14 Décembre 1991, puis du 14 Janvier au 17 Février permettant 30 vols à partir de l'aéroport de Kiruna. Le résultat le plus important obtenu est l'observation de la formation de nuages stratosphériques orographiques de grande extension horizontale (600 km) correspondant à des taux de refroidissement adiabatique de 3 à 10 K induits par les ondes de relief. La présence répétée du vortex polaire près de la Scandinavie, a permis de constater, à l'instar des mesures lidar de Sodankyla, la faible pénétration du nuage volcanique à haute altitude dans le vortex [Godin et al., 1994b].
### 3.4. Données satellitaires

Le Centre National de Recherches Météorologiques de Météo-France recevait dans un délai de quelques heures les mesures du sondeur TOVS émises par le Service opérationnel Américain NOAA / NESDIS (Washington) sur un domaine couvrant l'Atlantique Nord et une bonne moitié du continent Eurasien (20° à 90°N, 75°O à 90°E). Deux champs journaliers d'ozone total (orbites descendantes à 8h locales et montantes à 18h locales à calculées avec l'algorithme du CNRM, étaient disponibles le lendemain à 7h sur la base de données NADIR. Malgré un défaut d'étalonnage de l'algorithme dû à un changement de satellite (NOAA-10 remplacé par NOAA-12 en Septembre 1991), la structure fine du champ d'ozone total a pu être restituée avec une résolution horizontale de 20 à 60 km durant toute la campagne. Un algorithme de détection des nuages semi-transparents à en outre permis de localiser les nuages stratosphériques.

### 3.5. Modèles

Trois groupes français ont conduit des activités de simulation nu- mérique pour interpréter les résultats d'observation.

3.5.1. *Un modèle photochimique de boite* le long de trajec- toires Lagrangienne a été utilisé par le Service d'Aéronomie, le Laboratoire de Physique et Chimie de l'Environnement et l'ONERA pour simuler les mesures de distribution verticale de NO2 obtenues à l'aide des sondes

SAOZ. Il a été montré que les très faibles concentrations de NO₂ rapportées au mois de Janvier durant la partie la plus froide et la plus faiblement éclairée de l'hiver étaient simplement dues à l'absence de lumière solaire et pouvaient être simulées convena- blement en ne tenant compte que des réactions chimiques en phase gazeuse. Par contre, à température plus élevée en automne et au printemps, la simulation des faibles concentrations observées requière l'addition d'une conversion des NOx en acide nitrique. La réaction responsable est la conversion de N₂O₅ en HNO₃ sur la surface des gouttelettes d'aérosols volcaniques sulfatés; la réaction qui implique ClONO₂ est négligeable. Cependant la conversion n'est pas totale et les résultats montrent qu'une saturation du mécanisme hétérogène est observée qui provient de la formation très lente de N₂O₅ à partir de NO₂ à basse température [Lateltin et al., 1994].

3.5.2. *Le modèle de circulation générale MCG ARPEGE du Centre National de Recherches Météorologiques de Météo-France* a été utilisé pour étudier l'influence des condition météorologiques de l'hiver 1991-92 sur l'évolution de la couche d'ozone.

Un schéma d'initialisation du champ d'ozone à partir d'une combinaison de données satellitaires et de profils verticaux mesurés par les nombreux sondages, a été employé pour simuler les "mini trous" d'ozone qui se sont présentés au dessus de la Scandinavie vers la fin du mois de Janvier 1992. Il a été ainsi démontré que l'apparition de régions de faible colonne d'ozone à haute latitude telles que celles observées durant EASOE, étaient essentiellement le résultat d'une advection d'air pauvre en ozone en provenance des latitudes subtropicales [Orsolini et al., 1994].

3.5.3. *Un modèle de chimie-transport tridimensionnel forcé* par des champs de vent et de température provenant du Centre Européen de Prévisions Météorologiques à Moyen Terme (CEPMT), a également été utilisé pour évaluer l'étendue de l'action photochimique des nuages stratosphériques (PSC) et le potentiel de destruction d'ozone qui lui était lié, au cours de l'hiver 1991-92. Ceci est obtenu en calculant d'une part la durée durant laquelle l'air est exposé aux PSC, puis celle durant laquelle le même air est exposé au soleil. Il a été montré que l'action des PSC est d'abord demeurée épisodique au mois de décembre, puis s'est étendue à l'ensemble de l'Arctique au début

du mois de Janvier, pour s'évanouir par la suite lors des premiers réchauffements vers la fin du mois de Janvier. La conséquence de ceci est que les pertes d'ozone, bien que sensibles à certaines des approximations du modèle, sont demeurées faibles durant EASOE, à cause de l'absence de soleil durant la période froide où les PSC étaient présents [Chipperfield et al., 1994a].

L'étendue du mécanisme d'activation des chlores et des bromes a également été étudiée à l'aide du même modèle. La conclusion est que si la concentration de chlore actif a été très élevée (1.5 ppbV) durant la période froide au cours de la première partie du mois de Janvier, principalement à cause de la formation de PSC (et non à cause des aérosols volcaniques), elle a décru rapidement dès l'apparition des premiers échauffements [Chipperfield et al., 1994b].

3.5.4. *Le Laboratoire de Météorologie Dynamique* a cherché à interpréter les structures fines des profils verticaux d'ozone et de vapeur d'eau obtenus durant EASOE. Deux cas ont été étudiés, en particulier le jour où le lidar embarqué sur l'ARAT à mis en évidence l'existence d'un nuage stratosphérique qui portait la trace d'une onde de relief [Godin et al., 1994b]. C'est ainsi qu'il a été montré au moins dans les deux cas analysés, que les structures à petite échelle de la concentration des espèces, étaient la conséquence d'un transport induit par la propagation d'ondes de gravité [Teitelbaum et al., 1994].

## 4. CONCLUSION

Les résultats de l'ensemble de la campagne ont été publiés dans une édition spéciale du Geophysical Research Letters parue au mois de Juin 1994, à laquelle on se reportera pour plus de détails. Les principaux enseignements retirés de EASOE peuvent se résumer ainsi.

-La stratosphère de l'hémisphère Nord a été profondément perturbée en hiver 1991-92. A une atmosphère déjà riche en chlores et bromes d'origine industrielle, est venu s'ajouter une grande quantité d'aérosols volcaniques injectés à haute altitude lors de l'éruption du Pinatubo au mois de Juin précédent. Ceci a conduit à une situation nouvelle, jamais observée auparavant, qui est interprétée comme la cause de l'affaiblissement sensible des colonnes d'ozone observées par la suite en 1992 et 1993.

-Les aérosols ont été observés de façon extensive durant EASOE par plusieurs techniques qui toutes indiquent qu'il n'avaient pas ou peu pénétré dans le vortex polaire qui s'était formé à l'automne 1991.

-Les mesures de composés chlorés, ClO, OCIO, ClONO₂ et HCl, y compris les premières observations de la distribution verticale de OCIO et de ClONO₂ par des instruments embarqués en ballon, indiquent que les chlores se trouvaient sous formes actives durant la période froide au mois de Janvier dans le vortex polaire, mais que la concentration de celles-ci a rapidement diminuée au mois de Février lors du réchauffement de la stratosphère.

-Les mesures de distribution verticale et de colonne de NO₂ montrent que les NOx étaient quasiment absents à l'altitude de aérosols volcaniques. Ils étaient convertis en acide nitrique par réaction hétérogène sur les aérosols.

-Les quantités d'ozone étaient anormalement basses durant la campagne : -10% à -15% par rapport à la moyenne climatique de Décembre à Février. Cependant cet affaiblissement n'est pas la conséquence d'une destruction photochimique. Il est expliqué par la situation météorologique inhabituelle qui a prévalu en Europe au cours de l'hiver 1991-92 (hautes pressions record depuis 100 ans).

-Les taux de destruction d'ozone simulés par les modèles qui reproduisent bien à la fois l'effondrement des NOx et l'activation des chlores observés, ont été parfois importants durant EASOE.
Leur étendue a cependant été limitée dans l'espace et dans le temps et le phénomène s'est rapidement arrêté avec le réchauffement de la stratosphère au début du mois de Février.

En conclusion, les données de EASOE ont montré que le potentiel de destruction d'ozone par les chlores en Arctique était important. Du fait de la date du réchauffement stratosphérique, la destruction photochimique de l'ozone a été faible en hiver 1991-92, mais il apparaissait vraisemblable qu'avec l'inévitable augmentation de la charge en chlore dans les dix ans à venir, on pourrait assister à une destruction importante de l'ozone en Arctique en cas d'hiver froid et prolongé.

## REFERENCES

Aimedieu, P. Y. Kondo, T. Sugita, M. Koike and R. Sheldon, In situ measurements of nitric oxide, ozone and aerosol in the Scandinavian arctic stratosphere in January 1992, *Geophys. Res. Lett.*, **13**, 1243, 1994.

Chipperfield, M., D. Cariolle and P. Simon, A three-dimensional transport model study of PSC processing during EASOE, *Geophys. Res. Lett.*, **13**, 1463, 1994a.

Chipperfield, M., D. Cariolle and P. Simon, A 3D transport model study of chlorine activation during EASOE, *Geophys. Res. Lett.*, **13**, 1467, 1994b.

Flesia, C., A. Mugnai, Y. Emery, S. Godin, L. de Schoulepnikoff and V. Mitev, Interpretation of lidar depolarization measurements of the Pinatubo aerosols inferred with a microphysical code : A case study during EASOE, *Geophys. Res. Lett.*, **13**, 1435, 1994.

Godin, S., G. Mégie, C. David, D. Haner, Y. Emery, C. Flesia, V. Rizi, G. Visconti, L. Stefanutti, Proc. Quad. Ozone Symp., NASA Conf. Pub. 3266, 1994a.

Godin, S., G. Mégie, C. David, D. Haner, C. Flesia and Y. Emery, Airborne lidar observation of mountain-wave induced polar stratospheric clouds during EASOE, *Geophys. Res. Lett.*, **13**, 1335, 1994b.

Goutail F., J.P. Pommereau, A. Sarkissian, E. Kyro and V. Dorohkov, Total nitrogen dioxide at the Arctic polar circle since 1990, *Geophys. Res. Lett.*, **13**, 1371, 1994.

Herman, M., C. Brogniez, P. Lecomte, C. Verwaerde, Aerosol and PSC observations in polarized light during the EASOE campaign, *Proc. Second European Symp. on polar stratospheric ozone*, 1992.

Kondo, Y., U. Schmidt, T. Sugita, P. Aimedieu, M. Koike, H. Ziereis and Y. Iwasaka, Total reactive nitrogen, N2O and ozone in the winter Arctic stratosphere, *Geophys. Res. Lett.*, **13**, 1247, 1994.

Lateltin E., J. P. Pommereau, H. LeTexier, M. Pirre, R. Ramaroson, Perturbation of stratospheric nitrogen dioxide by volcanic aerosol in the arctic, *Geophys. Res. Lett.*, **13**, 1411, 1994.

Murcray, F.J., J.R. Starkey, W.J. Williams, W.A. Matthews, U. Schmidt, P. Aimedieu and C. Camy-Peyret, HNO3 profiles obtained during the EASOE Campaign, *Geophys. Res. Lett.*, **13**, 1223, 1994.

Neuber, R., G. Beyerle, G. Fiocco, A. di Sarra, K.H. Fricke, Ch. David, S. Godin, B.M. Knudsen, L. Stefanutti, G. Vaughan and J-P. Wolf, Latitudinal distribution of stratospheric aerosols during the EASOE winter 1991/92, *Geophys. Res. Lett.*, **13**, 1283, 1994.

Orsolini, Y., D. Cariolle and M. Déqué, A GCM study of the late January 1992 "mini-hole" event observed during EASOE, *Geophys. Res. Lett.*, **13**, 1459, 1994.

Ovarlez, J. and H. Ovarlez, Stratospheric water vapor content evolution during EASOE, *Geophys. Res. Lett.*, **13**, 1235, 1994.

Pommereau, J.P., Preliminary results of the European Arctic Ozone Experiment, *NATO Adv. Study Institute, Springer Ed.*, 1992.

Pommereau, J. P., H.K. Roscoe, J.A. Pyle, G.O.Braathen, S. Joffre, T.S. Jørgensen, M.Gil, P.C. Simon, G. Giovanelli and D. Perner, Interpretation of Stratosphere Monitoring by Ground-based Ultraviolet and Visible Spectrometry, *Rapport final à la Commission Européenne, Contrat STEP 0013M*, 1993.

Pommereau, J. P., C. Camy-Peyret, J. Ovarlez, G. Fiocco, J. de La Noé, P.C. Simon, L.Delbouille, T.S., Jørgensen, M. Gil, G. Giovanelli, G.O. Braathen, S. Joffre, C. Zerefos, L. Stefanutti, C. Varostos, L. Wöste, J.P. Wolf, K. Henriksen, European Arctic Stratospheric Ozone Experiment, Ground-Based Observations and Related Measurements, STEP 0141 MNLA, *Rapport final à la Commission Européenne*, 1993.

Pommereau, J. P and J. Piquard, Stratospheric Minor Species vertical Distributions During Polar Winter by Balloon borne UV-Visible Spectrometry, *Proc. Quad Ozone Symp.*, 528, NASA Conf. Pub. 3266, 1994.

Pommereau, J.P., F. Goutail, J. Piquard, L. Denis et C. Phillips, The SAOZ balloon optical sonde for atmospheric chemistry study, *Proc. 11th ESA Symp. Eur. Rocket & Balloon Prog.*, ESA SP 355, 87, 1994.

Pommereau, J. P. and J. Piquard, Ozone, Nitrogen dioxide and Aerosol vertical distributions by UV-visible solar occultation from balloons, *Geophys. Res. Lett*, **13**, 1227, 1994.

Pommereau, J.P., and J. Piquard, Observations of the vertical distribution of stratospheric OCIO, *Geophys. Res. Lett.*, **13**, 1231, 1994.

Pyle, J.A, N.R.P. Harris, J.C. Farman, F. Arnold, G. Braathen, R.A. Cox, P. Faucon, R.L. Jones, G. Megie, A.O'Niel, U. Platt, J.P. Pommereau, U. Schmidt and F. Stordal, An Overview of the EASOE Campaign, *Geophys. Res. Lett.*, **13**, 1191, 1994.

Sarkissian A., J. P. Pommereau, F. Goutail and E. Kyro, PSC and volcanic aerosol observations during EASOE by UV-visible ground-based spectrometry, *Geophys. Res. Lett.*, **13**, 1319, 1994.

Teitelbaum, H., J. Ovarlez, H. Kelder and F. Lott, Some observations of gravity-wave-induced structure in ozone and water vapour during EASOE, *Geophys. Res. Lett.*, **13**, 1483, 1994.
# IMAGERIE ELECTROMAGNETIQUE BAYESIENNE

Michel MENVIELLE

Équipe de Physique de la Terre et des Planètes (URA 1369), Bâtiment 504, Université Paris Sud, F-91405 ORSAY Cedex

Michel ROUSSIGNOL

Équipe d'Analyse et de Mathématiques Appliquées, Université de Marne la Vallée, 2 rue de la Butte Verte, 93166 NOISY LE GRAND

Cet article présente la contribution des équipes françaises au développement de méthodes d'imagerie électromagnétique Bayesienne. Les méthodes développées par ces équipes sont basées sur deux approches différentes, l'agorithme de Monte Carlo d'une part, et l'utilisation d'une chaîne de Markov ergodique dont la loi de transition se stabilise asymptotiquement vers la loi a posteriori des paramètres. Ces méthodes permettent d'échantillonner de façon satisfaisante cette loi a posteriori pour une une grande variété de lois a priori, et ce sans qu'il soit nécessaire de linéariser la fonction du problème direct.

This paper presents the results established by french teams in Bayesian electromagnetic imaging. The methods they developed are based on the one hand on Monte Carlo algorithm, and on the other hand on the use of an ergodic Markov chain the transition law of which tends asymptotically towards the a posteriori law of the parameters. These methods allow to estimate the a posteriori law for a wide class of a priori laws. The computation are made without linearization of the forward problem function.

## 1. INTRODUCTION

Notre compréhension actuelle de la structure et de la dynamique interne de la Terre résulte de travaux pluridisciplinaires qui combinent les résultats de nombreuses disciplines des sciences de la Terre, telles que la pétrologie, la géochimie, la géophysique, l'étude en laboratoire des propriétés physiques et chimiques des matériaux terrestres, la simulation numérique... Les résultats fournis par chacune de ces disciplines constituent en effet autant d'informations indépendantes, qui se complètent et se précisent mutuellement, permettant d'estimer les conditions qui règnent à l'intérieur de la planète.

Une des limitations majeure de notre compréhension de la structure et la dynamique de la Terre provient ainsi de la précision avec laquelle sont connues les distributions de paramètres physico-chimiques indépendants, tels que la vitesse de propagation des ondes élastiques, la densité, de la susceptibilité magnétique et la résistivité électrique pour ne citer que certains de ces paramètres. Ces distributions sont déduites de l'interprétation de mesures géophysiques effectuées en surface, ou à bords d'aéroneks et de satellites équipés de capteurs appropriés. Les techniques de mesure et de traitement des données ont fait des progrès remarquables ces dernières années, et elles permettent maintenant de disposer de données d'excellente qualité. La précision des méthodes quantitative d'interprétation des observations est ainsi ce qui limite le plus souvent notre connaissance de la Terre, et le développement de méthodes précises et efficaces d'imagerie géophysique est un des enjeux majeurs des années à venir.

C'est en particulier le cas en électromagnétisme, où ces méthodes ne sont réellement arrivés à maturité que pour des situations présentant des particularités géométriques remarquables (milieux tabulaires où la conductivité ne dépend que de la profondeur, ou structures d'extension infinie dans une des directions horizontales). L'imagerie de structures tridimensionnelles, fondamentale à bien des égards, est encore balbutiante. L'amélioration des méthodes d'imagerie électromagnétique est indispensable pour progresser dans notre connaissance de la structure et de la dynamique interne de la Terre, car la conductivité des matériaux terrestres est très étroitement reliée à leur pétrologie et aux conditions thermodynamiques qui prévalent dans le milieu étudié [Olhoeft, 1980; Parkhomenko, 1982].

Le champ électromagnétique est régi par les équations de Maxwell, et le problème de l'imagerie électromagnétique est un problème non linéaire, qui n'admet une solution unique que pour des situations jamais rencontrées en pratique en géophysique. Résoudre le problème inverse en électromagnétisme nécessite donc de caractériser l'ensemble des modèles acceptables, qui permettent d'expliquer les observations, et si possible quantifier les informations qu'apportent les observables étudiées sur notre connaissance de ces modèles (voir par exemple [Backus, 1988] pour une discussion de la non unicité de la solution).

L'approche Bayesienne constitue un cadre très adapté à la description de l'ensemble des modèles acceptables, car elle associe à chaque modèle une probabilité a posteriori qui décrit sous forme probabiliste la compatibilité de ce modèle avec les observations et l'information sur la zone étudiée éventuellement disponible par ailleurs. Cet article présente la contribution des équipes françaises au développement de méthodes d'imagerie électromagnétique Bayesienne. Nous rappelons d'abord brièvement les diverses méthodes de résolution du problème inverse non linéaire, puis nous développons plus particulièrement les méthodes Bayesienne et présentons les résultats établis dans ce domaine par des équipes françaises.

## 2. LE PROBLEME INVERSE NON-LINEAIRE

La résolution du problème inverse a deux objectifs :

- trouver un modèle acceptable, qui permette de rendre compte des observables, et, si possible, caractériser l'ensemble des modèles acceptables, d'une part;

- quantifier l'information apportée par les données traitées, d'autre part.

Ces objectifs sont en fait très difficiles à atteindre en totalité lorsque le problème est non linéaire, comme c'est le cas en électromagnétisme, et l'objectif de nombreuses méthodes de résolution du problème inverse se réduit en fait à trouver le modèle qui rende compte au mieux des observations, ou observables.

Le principe de base de la résolution du problème inverse consiste à chercher un modèle dont la réponse théorique est la plus proche possible des observables. En d'autres termes, la différence entre les données observées et celles calculées à l'aide d'un modèle est estimée à l'aide d'une fonction d'écart, et le
modèle recherché est celui qui minimise cette fonction. D'un point de vue topologique, le choix de la fonction d'écart définit une distance dans l'espace des modèles, et la validité d'un modèle est mesurée par la distance entre ce modèle et les observables; d'un point de vue probabiliste, la fonction d'écart définit la probabilité d'un modèle, probabilité qui est d'autant plus grande que la fonction d'écart est petite.

Les méthodes de résolution du problème inverse non-linéaire peuvent être divisées selon la manière dont on cherche la solution dans l'espace des modèles :

- une première famille est basée sur une approche globale. Un sous ensemble de l'espace des modèles, correspondant à des valeurs des paramètres comprises dans des intervalles bornés, est déterminé à partir du savoir *a priori* sur les valeurs possibles des paramètres du modèle. Ce sous ensemble est exploré à l'aide d'un processus aléatoire permettant de générer une série de modèles décrivant ce sous ensemble avec une précision donnée, et la solution du problème inverse est déterminée à partir des modèles ainsi obtenus. Cette approche, généralement appelée méthode de Monte Carlo [Metropolis et Ulam, 1949; Andersen et Seneta, 1971; Tarantola, 1987], permet d'éviter la linéarisation de la fonction d'écart, et les problèmes liés à cette approximation. Mais elle nécessite en contrepartie une puissance de calcul considérable, car l'exploration d'un sous ensemble de l'espace des modèles correspondant à des hypothèses réalistes sur les paramètres impose de générer un grand nombre de modèles;

- une seconde famille de méthodes a ainsi été développée afin d'obtenir une solution acceptable avec un coût raisonnable en temps calcul. Ces méthodes sont basées sur une approche locale. Etant donné un modèle de départ, choisi *a priori* et traduisant en général la connaissance dont on dispose avant inversion sur le milieu étudié, le processus chemine dans l'espace des modèles par perturbations successives du modèle de façon à minimiser la fonction d'écart. En chaque point du cheminement, la direction à prendre est déterminée à partir du comportement de la fonction d'écart en ce point. Dans une approche déterministe, la direction choisie est celle qui optimise la décroissance de la fonction d'écart : c'est par exemple le cas des méthodes de gradient maximum, dans lesquelles la direction choisie est celle qui optimise la décroissance en module de la fonction d'écart linéarisée. Dans une approche probabiliste, le choix de la direction est effectué à partir d'une distribution de probabilité sur les directions qui dépend du comportement de la fonction d'écart autour du point : c'est par exemple le cas des méthodes de recuit simulé (cf. par exemple [Dosso et Oldenburg, 1991], ou de celles basées sur l'algorithme de Metropolis (cf. par exemple [Metropolis et al., 1953; Marroquin et al., 1987; Mosegaard et Tarantola, 1995]). Dans tous les cas, la solution est le modèle qui correspond au minimum de la fonction d'écart dans laquelle se stabilise le processus (cf. par exemple [Wu, 1968; Pedersen et Hermance, 1986; Tarantola, 1987; Pedersen et Rasmussen, 1989; Fisher et Le Quang, 1981; Menke, 1984; Press et al., 1989]; mais il faut cependant noter que les méthodes de cheminement basées sur une approche déterministe (méthodes de gradient maximal par exemple) ne permettent pas d'empêcher le piégeage dans un minimum secondaire de la fonction d'écart, alors que les méthodes basées sur une marche au hazard (Monte Carlo, Metropolis, chaîne de Markov, pour ne citer que les plus connues) permettent d'échantilloner de façon satisfaisante la fonction d'écart si un nombre suffisant d'itérations est effectué. Enfin, de nouvelles méthodes stochastiques simulant des processus existant dans la nature, tels que l'algorithme génétique (cf. par exemple [Sambridge et Drijcknningen, 1992; Everett et Schultz, 1993], ont été proposées ces dernières années;

- une troisième famille de méthodes, enfin, est basée sur

l'emploi d'opérateurs mathématiques, en général des transformations intégrales. Ces opérateurs permettent de calculer, à partir du champ électromagnétique observé, des grandeurs plus aisément interprétables qui seront utilisées pour caractériser la structure conductrice du milieu sous jacent (voir [Spichak, 1990] pour une présentation synthétique de ces méthodes). Une classe particulière de transformations intégrales du champ électromagnétique correspond aux méthodes asymptotiques, développées par analogie avec les méthodes de rétropropagation des ondes sismiques. Ces méthodes sont basées sur des opérateurs permettant d'associer un champ électromagnétique, solution de l'équation hyperbolique de la diffusion, à une fonction solution de l'équation parabolique de la propagation des ondes.

La contribution des équipes françaises à la résolution du problème inverse en électromagnétisme concerne chacun des aspects du problème inverse que nous venons d'évoquer :

- le développement de méthodes Bayesiennes permettant de déterminer le modèle le plus probable, et de décrire de façon probabiliste l'information disponible après inversion [Tarits et Jouanne, 1990; Jouanne, 1991; Roussignol et al., 1993; Grandis, 1994; Tarits et al., 1994; Grandis et al., 1995a; Grandis et al., 1995b; Menvielle et Roussignol, 1995; Spichak et al., 1995];

- le développement de méthodes basées sur l'utilisation de méthodes asymptotiques [Gibert et Virieux, 1991; Gibert et al., 1994; Virieux et al., 1994; Tournerie et al., 1995].

Le présent article a pour objectif la présentation des résultats concernant le développement de méthodes Bayesiennes de résolution du problème inverse en électromagnétisme.

### 3. IMAGERIE ELECTROMAGNETIQUE BAYESIENNE

Dans une perspective statistique Bayesienne, la résolution du problème inverse consiste à supposer que l'information disponible sur le modèle après inversion, le savoir *a posteriori*, résulte de l'actualisation de l'information disponible avant inversion, le savoir *a priori*, à l'aide de toute celle apportée par les observables. Les savoirs *a priori* et *a posteriori* sont exprimés sous forme probabiliste, par les lois *a priori* et *a posteriori* des paramètres. Résoudre le problème inverse, c'est alors déterminer la loi *a posteriori* des paramètres sachant leur loi *a priori* et celle des observables, c'est à dire exprimer sous forme probabiliste l'actualisation, à l'aide des informations contenues dans les données interprétées, du savoir *a priori* sur le modèle [Box et Tiao, 1973; Berger, 1985; Press, 1989].

La résolution du problème inverse nécessite souvent de faire un certain nombre d'hypothèses, explicites ou non, sur le modèle recherché : caractéristiques géométriques du milieu étudié, domaine de variation des paramètres, caractère plus ou moins lisse des variations spatiales des paramètres, pour ne citer que les plus fréquentes. Parallèlement à ces hypothèses, il est souvent nécessaire d'effectuer un certain nombre de choix qui conditionnent le résultat sans qu'il soit toujours aisé d'en quantifier l'effet. Ainsi, dans le cas des méthodes de cheminement déterministe, l'influence sur le résultat du choix du modèle initial n'est jamais aisée à déterminer. L'approche Bayesienne fournit un cadre qui permet d'explicitier les hypothèses faites, au moyen de la loi *a priori* des paramètres, et de quantifier l'information apportée dans ce contexte par les observables, en comparant la loi *a posteriori* avec la loi *a priori*.

La résolution du problème inverse n'est possible que si le problème direct est résolu. Nous noterons dans la suite F la fonction du problème direct, qui permet de calculer la valeur de l'observable y lorsque le modèle est x. En considérant que l'erreur e provient uniquement de l'observation, les données
actuellement observées sont :

$$
y = F(x) + \epsilon \tag{3.1}
$$

La probabilité *a posteriori* des paramètres **x** sachant les observables *y* est exprimée par la formule de Bayes [Bayes, 1763] :

$$
p(\mathbf{x}/\mathbf{y}) = \frac{p(\mathbf{y}/\mathbf{x}) g(\mathbf{x})}{\int_{\mathbf{x}} p(\mathbf{y}/\mathbf{x}) g(\mathbf{x}) \, \mathrm{d}\mathbf{x}} \quad (3.2)
$$

où p(a/b) est la probabilité conditionnelle de a sachant b, et
g(x) la loi a priori des paramètres x. Dans le cas général, il est en
pratique impossible de calculer directement la loi a posteriori
p(x/y), car l'expression de cette loi comporte des intégrales
multiples dont le noyau est une fonction compliquée des
paramètres.

L'expression explicite de cette loi n'est en fait possible que
lorsque les lois sur l'observable et sur les paramètres sont
gaussiennes, et que la fonction F est linéaire en fonction des
paramètres du modèle [Menke, 1984; Tarantola, 1987; Backus,
1988]. Dans le cas d'un problème inverse non linéaire, comme
c'est le cas pour l'électromagnétisme, cette possibilité n'est pas
très intéressante en pratique car elle impose de linéariser F, avec
tous les risques d'instabilités que cela comporte. De plus
l'hypothèse d'une loi gaussienne sur les paramètres implique que
la valeur la plus probable de ces paramètres soit connue. (voir
par exemple [Mackie et al., 1988]. Cette hypothèse est contrai-
gnannte en électromagnétisme, car le paramètre étudié, la
résistivité des matériaux géologiques, peut varier de plusieurs
ordres de grandeur entre quelques milliers d'Ωm (roche froide
sèche) et quelques dixièmes d'Ωm (roche contenant un fluide
conducteur : fluide hydraté ou roche fondue).

3.1. Méthodes de Monte Carlo

Les méthodes d'intégration stochastiques du type Monte Carlo
constituent, dans un certain nombre de cas, un outil relativement
bien adapté à l'estimation de la loi *a posteriori* des paramètres
donnée par la formule (3-2).

Une solution pour optimiser le calcul des intégrales multiples par la méthode de Monte Carlo est de sélectionner à l'aide d'une fonction appelée fonction d'importance les valeurs pour lesquelles l'intégrale est estimée. Cette fonction est en pratique une densité de probabilité choisie de manière à favoriser des échantillons dont la contribution à l'intégrale est la plus importante. Utilisant une méthode proposée par [Wiggins, 1969; Wiggins, 1972], [Jouanne, 1991] et [Tarits et al., 1994] ont utilisé la loi a priori comme fonction d'importance pour le calcul de la loi a posteriori. Ces auteurs ont ainsi développé une méthode de Monte Carlo qui permet en pratique d'estimer numériquement les lois marginales a posteriori des paramètres pour un problème inverse non linéaire, dès l'instant où le nombre de paramètres du modèle reste petit (typiquement 10-20 paramètres), et ce sans qu'il soit nécessaire de faire d'hypothèses contraignantes sur la loi a priori des paramètres.

[Jouanne, 1991] et [Tarits et al., 1994] ont utilisé cette
approche pour résoudre le problème inverse magnétotellurique
pour un milieu sous jacent tabulaire (situation 1-D) composé de
quelques couches homogènes au dessus d'un demi espace
homgène lui aussi. Les paramètres du modèle sont les résistivités
de chacun des milieux homogènes et les profondeurs des limites
entre ces milieux. Le savoir a priori concerne donc le nombre de
couches, et le domaine dans lequel chacun des paramètres peut
varier, avec une loi a priori uniforme.

Ces auteurs ont utilisé ce modèle pour discuter l'influence

d'erreurs dans le savoir *a priori*. Une erreur dans la détermination du domaine de variation possible d'un paramètre est en général très facile à détecter, car elle se traduit par l'existence d'un maximum de la loi marginale *a posteriori* correspondante sur l'une des bornes du domaine de variation. Ils ont aussi étudié l'influence d'une erreur sur la détermination du nombre de couche. Un nombre de couches trop élevé se traduit dans la morphologie des lois *a posteriori*, et l'analyse de ces lois permet de déterminer le nombre de couches adapté à l'information contenu dans les données. En revanche, les lois *a posteriori* obtenues pour un nombre de couches trop faible ne permettent pas de détecter une telle erreur, et le résultat obtenu peut être significativement biaisé. Ce résultat met clairement en évidence le rôle crucial de la détermination *a priori* du nombre de couches, à partir, par exemple, du résultat d'inversions préliminaires à l'aide de méthodes permettant d'estimer simplement ce nombre, comme celle du D⁺ proposée par [Parker, 1980].

[Tarits et Jouanne, 1990] ont utilisé cette approche pour étudier quantitativement la variation de la structure de la lithosphère océanique avec l'âge à partir de l'analyse de sondages électromagnétiques fond de mer effectués sur des sites d'age différent dans le Pacifique Nord [Filloux, 1977; Filloux, 1980; Law et Greenhouse, 1981; Yukutake et al., 1983; Tarits. 1989]. Les profils de variation de la résistivité électrique en fonction de la profondeur obtenus par [Tarits et Jouanne, 1990] sont présentés dans la figure 1a. Ils mettent en évidence une dépendance linéaire entre la profondeur de la limite entre une "lithosphère électrique" résistante et une "asthénosphère électrique" conductrice d'une part, et la racine carrée de l'âge de la lithosphère d'autre part (Figure 1b)

3.2. Méthodes Markoviennes

Du fait de leur coût élevé en temps calcul, les méthodes de
Monte Carlo deviennent rapidement impossible à mettre en
oeuvre dès que le nombre de paramètres augmente.

[Jouanne, 1991; Roussignol et al., 1993; Grandis, 1994] ont développé une méthode de résolution du problème inverse basée sur l'utilisation d'une chaîne de Markov ergodique dont la loi de transition se stabilise asymptotiquement vers la loi a posteriori des paramètres. Ils ont appliqué cette méthode à la résolution du problème inverse électromagnétique pour des structures de géométrie variée : tabulaire (cas unidimensionnel), planes et d'extention latérale grande devant leur épaisseur (plaques minces), allongée dans une direction (cas bidimensionnel). Après avoir exposé les bases théoriques de la méthode de résolution du problème inverse proposée par ces auteurs, nous présenterons les principaux résultats obtenus par ces auteurs pour des données sythétiques et réelles.

3.2.1. Théorie. Nous allons, dans un premier temps, rappeler les propriétés fondamentales des chaînes de Markov, puis nous introduirons les algorithmes proposés par ces auteurs en suivant la démarche adoptée par [Roussignol et al., 1993], et enfin nous présenterons les principaux résultats théoriques établis par ces auteurs.

Considérons un système qui peut prendre un certain nombre
d'états, et dont l'évolution au cours du temps correspond à un
processus aléatoire. A un instant donné, l'état du système peut
donc être décrit comme une variable aléatoire appartenant à
l'ensemble E des états possibles. L'évolution du système est un
processus Markovien si, à chaque instant, l'évolution future du
système ne dépend que de son état présent, ce qui signifie que le
système ne dépend du passé qu'au travers de son état actuel. Les
chaînes de Markov constituent la classe particulière de processus
de Markov tels que l'espace des états possibles est fini ou
dénombrable, et que les instants successifs du processus de
<figure><img src="image_110.png" /></figure>

Figure 1.a. Profils de résistivité obtenu par inversion à l'aide de la méthode développée par [Tarits et al., 1994] de sondages fond de mer effectués sur des sites d'age différents dans le Pacifique Nord : JDF : 1Ma; EMS : 5Ma; CAL : 35 Ma; HAW : 70 Ma; NCP : 72 Ma; STJ : 130 Ma. (d'après [Tarits et Jouanne, 1990]).
<figure><img src="image_111.png" /><figcaption>Figure 1.b. Variation de la profondeur du toit de l'asthénosphère conductrice en fonction de la racine carrée de l'age de la lithosphère (d'après [Tarits et Jouanne, 1990])</figcaption></figure>

Markov sont repérés par des entiers. Une suite {$X_i$; i=0, 1, 2, ... n} de variables aléatoires à valeurs dans un espace fini ou dénombrable E est donc une chaîne de Markov si l'état du système au temps n+1 dépend du passé uniquement par l'intermédiaire de l'état du système au temps n. Le comportement de la chaîne de Markov {$X_i$} est défini par l'ensemble de ses probabilités de transition :

$$P_n(a, b) = P(X_{n+1} = b | X_n = a), a, b \in E \quad (3.3)$$

Lorsque (i) les probabilités de transition ne dépendent pas du temps n, (ii) E est fini, et (iii) tous les états possibles sont accessibles chacun de l'un à l'autre, la chaîne est dite ergodique. Lorsque n augmente indéfiniment, le comportement asymptotique des chaînes de Markov ergodiques est tel que :

-la moyenne empirique des lois de $X$ converge vers la loi invariante $\Pi$ lorsque n tend vers l'infini, et ce indépendamment de l'état initial. Cette moyenne est donc un estimateur de la loi de probabilité invariante $\Pi$;

-le nombre moyen de fois où la chaîne est dans un état $b$ converge vers la probabilité invariante de $b$, $\Pi(b)$ : la moyenne empirique des valeurs successives $X_n$ de la chaîne converge, lorsque le temps n tend vers l'infini, vers la moyenne probabiliste des éléments de E muni de la loi de probabilité $\Pi$.

Supposons maintenant que le problème inverse Bayesien consiste à déterminer les paramètres d'un modèle constitué de L éléments {$x_i$, i=1, 2, ..., L}, caractérisés chacun par un paramètre. Supposons en outre qu'il n'y ait que $M_l$ valeurs possibles de paramètre $\theta_{ml}$; m=1, 2, 3, ... $M_l$ associées à chaque élément $x_i$. Le choix des valeurs possibles $\theta_{ml}$ exprime le savoir a priori sur les paramètres. La probabilité marginale a posteriori que l'élément $x_i$ prenne la valeur $\theta_{ml}$ quelles que soient les valeurs des autres éléments s'écrit alors :

$$\Pi(x; x_1 = \theta_{ml} | y) = \frac{\sum_{x \in X_{l,m1}} f(y/x; x_1 = \theta_{ml}) g(x; x_1 = \theta_{ml})}{\sum_{x \in X} f(y/x) g(x)} \quad (3.4)$$

où X est l'ensemble des images possibles, et où {$x$; $x_i=\theta_{ml}$} et $X_{l,m1}$ représentent respectivement une image et l'ensemble des

images où $x_i=\theta_{ml}$. La constante de normalisation, qui apparaît au dénominateur de l'expression de la probabilité marginale a posteriori, est difficile, en pratique, à estimer car c'est une somme sur l'ensemble des images possibles, de dimension $M_1 \cdot M_2 \cdot M_3 \cdot ... \cdot M_l$.

En revanche, étant donnée une image quelconque dans l'ensemble $X_{l,m1}$, la probabilité conditionnelle a posteriori de $x_i=\theta_{ml}$ sachant les valeurs de tous les autres paramètres se calcule aisément :

$$P(x; x_1 = \theta_{ml}) = \frac{f(y/x; x_1 = \theta_{ml}) g(x; x_1 = \theta_{ml})}{\sum_{j=1}^{M_1} f(y/x; x_1 = \theta_{jl}) g(x; x_1 = \theta_{jl})} \quad (3.5)$$

La constante de normalisation de cette probabilité conditionnelle a posteriori est une somme sur l'ensemble des $M_1$ valeurs possibles de l'élément $x_i$ et son calcul devient alors tout à fait accessible pour un nombre raisonnable $M_1$ de valeurs possibles $\theta_{ml}$.

Les probabilités conditionnelles a posteriori exprimées par l'équation (3.5) définissent une loi de probabilité des paramètres pour l'élément $x_i$, et elle peuvent être utilisées dans un processus itératif pour actualiser la valeur de $x_i$ à l'aide d'un tirage au sort de loi $P(a; x_i = \theta_{ml})$ parmi les valeurs possibles $\theta_{ml}$. Lorsque l'élément à actualiser est choisi de façon aléatoire, c'est à dire par un tirage au sort uniforme parmi L éléments, la suite des images obtenues est une chaîne de Markov ergodique. La loi invariante vers laquelle converge la probabilité de transition est la probabilité a posteriori des paramètres sachant le savoir a priori et l'observable y [Roussignol et al., 1993; Grandis, 1994]. Si le choix de l'élément à actualiser est déterministe (cas du balayage régulier des éléments, par exemple), la chaîne de Markov constitué par la suite des modèles obtenues en actualisant successivement chaque élément n'est plus homogène, mais la suite des modèles $M_n$ obtenus à la fin de chaque cycle de balayage complet sur tous les éléments constitue toujours une chaîne de Markov ergodique qui possède la loi a posteriori des paramètres comme loi invariante [Roussignol et al., 1993; Grandis, 1994]. C'est cette dernière chaîne de Markov que ces auteurs ont utilisé pour proposer une solution du problème inverse en électromagnétisme, et nous nous intéresserons uniquement à cette chaîne dans la suite.

Le calcul de la probabilité de transition de la chaîne de Markov nécessite de résoudre complètement le problème direct $M_1$ fois pour chaque élément $x_i$, et donc $M_1+M_2+M_3+...+M_1$ fois entre deux modèles successifs $M_n$ qui constituent la chaîne de Markov. Le volume de calcul nécessaire pour résoudre le problème inverse est bien moins important que dans le cas de processus de Monte Carlo, et cette chaîne de Markov a été utilisée pour résoudre le problème inverse dans le cas de structures tabulaires. Ce volume de calcul reste cependant trop important dans le cas de structures de géométrie plus complexe, pour lesquelles la fonctionnelle F du problème direct (équation 3-1) ne s'exprime pas aussi simplement que dans le cas d'un milieu tabulaire. Afin de limiter le temps de calcul, le problème direct n'est pas résolu complètement à chaque étape du balayage, et le champ électromagnétique est calculé par la chaîne de Markov. La valeur des paramètres (c'est à dire la conductivité des domaines homogènes du modèle) et le champ électromagnétique sont alors estimés simultanément, et la chaîne de Markov modifiée comporte une variable (la conductivité) distribuée suivant une loi discrète et une autre variable (le champ électromagnétique) distribuée suivant une loi continue sur un ensemble de vecteurs réels. Il s'agit alors d'une chaîne de Markov à deux dimensions.

Les propriétés mathématiques, et notamment le comportement
asymptotique d'une telle chaîne de Markov n'avaient pas été étudiés lorsque les premiers résultats obtenus en électromagnétisme à l'aide de cette chaîne ont été présentés par [Jouanne et al., 1989]. [Jouanne, 1991], puis [Grandis, 1994] ont étudié empiriquement le comportement asymptotique de cette chaîne pour des situations représentatives de celles rencontrées en géophysique. [Grandis, 1994] a en particulier montré sur des exemples synthétiques que la loi de transition de la chaîne se stabilise asymptotiquement vers une loi invariante, qui correspond à la loi *a posteriori* des paramètres, si le champ électromagnétique est estimé, à chaque étape du balayage, avec une précision suffisante pour permettre une estimation raisonnable des probabilités de transition. L'étude mathématique du comportement asymptotique d'une chaîne de Markov à deux dimensions a été entreprise, à l'initiative de M. Roussignol. Cette chaîne de Markov, notée formellement ($\mathbf{X}_n, \Theta_n$), est à valeurs dans $Exl$, où $E$ est un espace fini et $I$ une partie de $\mathbb{R}^k$. Dans le cas où $E$ est un espace fini et $I$ un intervalle de $\mathbb{R}^+$, [Touijar, 1994] a ainsi démontré que :

- la suite $(\Theta_n)_{n \ge 0}$ est proche de la solution d'une équation différentielle dont l'attracteur est le champ électrique correspondant à la vraie distribution de conductivité dans le milieu sous jacent;

- l'algorithme converge sous certaines conditions.

Bien que le critère de convergence établi soit très sévère lorsque le nombre d'éléments de $E$ est grand, comme c'est le cas pour les images que nous considérons, ces résultats fondent théoriquement un résultat (l'existence de chaînes convergentes) qui avait été établi expérimentalement à partir d'expériences numériques et extrapolé à l'ensemble des situations géophysiques. L'interprétation physique de ces résultats suggère que la convergence de la chaîne est étroitement liée à celle du champ électrique vers sa valeur asymptotique, convergence dont on peut raisonnablement supposer qu'elle varie peu d'une situation géophysique à l'autre. Les résultats théoriques établis par D. Touijar sont une première étape pour justifier l'extrapolation à l'ensemble des situations géophysique des propriétés de convergence observées sur des exemples représentatifs.

### 3.2.2. Cas unidimensionnel.
Dans le cas unidimensionnel, [Grandis, 1994] a considéré un modèle constitué d'un grand nombre de couches minces d'épaisseur donnée *a priori* : les paramètres du modèle sont donc les résistivités de chacune des couches. Cette approche, qui revient en fait à échantilloner la courbe de variation de la résistivité en fonction de la profondeur, permet de s'affranchir du problème crucial, et sur lequel peu d'informations sont généralement disponibles, de la détermination *a priori* nombre de couche du modèle.

L'algorithme utilisé n'impose aucune hypothèse sur la loi *a priori* des paramètres : [Grandis, 1994] a utilisé une loi uniforme sur un intervalle borné. Cette loi est la même pour tous les paramètres, c'est à dire que le savoir *a priori* se résume à la détermination des valeurs minimales, $\rho_{min}$, et maximales, $\rho_{max}$, de la résistivité dans le milieu étudié. Les valeurs possibles de la résistivité, $\theta_m$, ne dépendent pas de la couche 1 et sont équivalantes en logarithme sur l'intervalle $\{\rho_{min}, \rho_{max}\}$. Les observables sont les valeurs de l'impédance magnétotellurique scalaire, chacune des observables (un nombre complexe $Z_1$ pour chaque période $T_1$) est indépendante, et le bruit $\epsilon$ est supposé gaussien.

[Grandis, 1994] a résolu le problème inverse unidimensionnel ainsi posé à l'aide de la chaîne de Markov que nous venons de décrire. Du fait de la relative simplicité des calculs, il est possible de résoudre complètement le problème direct, et la chaîne utilisée est une chaîne ergodique classique dont la loi invariante est la loi *a posteriori* des paramètres, les résistivités des couches. La solution obtenue est instable, du fait

des caractéristiques physiques inhérentes au problème direct pour le paramétrage choisi et du nombre de paramètres à estimer qui est grand devant le nombre des observables.

Afin de stabiliser la solution obtenue, [Grandis, 1994] a favorisé les modèles réguliers, pour lesquels la variation de résistivité entre deux couches consécutives est aussi faible que possible. Ce procédé de régularisation des modèles est classique (cf. par exemple [Tikhonov et Arsenin, 1977; Constable et al., 1987; Smith et Booker, 1988]). Dans le formalisme du problème inverse Bayesien, l'introduction d'une contrainte pénalisant l'irrégularité des modèles peut être effectuée en modifiant la loi *a priori*, l'effet de cette contrainte sur les paramètres étant décrit par la loi invariante de la chaîne dont les probabilités de transition sont construites en tenant compte de la loi *a priori* de régularisation. [Grandis, 1994] a montré que, étant donnée la matrice $M$ dont chacun des éléments $M(\theta_i, \theta_j)$ est la mesure de l'irrégularité entre deux couches successives de résistivité $\theta_i$ et $\theta_j$, ($M(\theta_i, \theta_j)>0$, avec $M(\theta_i, \theta_i)=1$), la chaîne de Markov dont les probabilités de transition sont :

$$ 
\begin{align*}
P(\mathbf{a}; \rho_1 = \rho_m) &= \frac{f(\mathbf{d}/\mathbf{a}; \rho_1 = \rho_m) w(\rho_m) M(\rho_m, \rho_2)}{\sum_{j=1}^{M} f(\mathbf{d}/\mathbf{a}; \rho_1 = \rho_j) w(\rho_j) M(\rho_j, \rho_2)} \\
P(\mathbf{a}; \rho_k = \rho_m) &= \frac{f(\mathbf{d}/\mathbf{a}; \rho_k = \rho_m) M(\rho_{k-1}, \rho_m) M(\rho_m, \rho_{k+1})}{\sum_{j=1}^{M} f(\mathbf{d}/\mathbf{a}; \rho_k = \rho_j) M(\rho_{k-1}, \rho_j) M(\rho_j, \rho_{k+1})} \\
P(\mathbf{a}; \rho_L = \rho_m) &= \frac{f(\mathbf{d}/\mathbf{a}; \rho_L = \rho_m) M(\rho_{L-1}, \rho_m) v(\rho_m)}{\sum_{j=1}^{M} f(\mathbf{d}/\mathbf{a}; \rho_L = \rho_j) M(\rho_{L-1}, \rho_j) v(\rho_j)} \tag{3-6}
\end{align*}
$$

a pour loi invariante :

$$ \Pi(\mathbf{a}) = \frac{v(\mathbf{a}) \bullet w(\mathbf{a})}{\sum_{\mathbf{a} = \rho_1, \dots, \rho_m} v(\mathbf{a}) \bullet w(\mathbf{a})} \quad (3-7) $$

où $v=\{v(\cdot)\}$ est un vecteur propre à droite et $w=\{w(\cdot)\}$ un vecteur propre à gauche associés à la plus grande valeur propre de $M$. [Grandis, 1994] a montré que, dans le cas qui nous intéresse ici, une mesure de l'irrégularité de la forme :

$$ M(\theta_i, \theta_j) = \exp\left(-\lambda \left(\log\left(\frac{\theta_i}{\theta_j}\right)\right)^2\right) \quad (3-8) $$

est la plus adaptée car la loi *a priori* est la même sur tout l'intervalle $\{\rho_{min}, \rho_{max}\}$ si l'échantillonnage est régulier en logarithme. La probabilité de transition s'écrit alors, pour le cas général où $1<l<l$ :

$$ P(\mathbf{a}; \rho_i = \theta_{ml}) = Cte \exp \left[ -\frac{1}{2} \sum_{i=1}^{2N} \left( \frac{d_i - F_i(\mathbf{a}; \rho_i = \theta_{ml})}{\sigma_i} \right)^2 - \lambda \sum_{i=k}^{k+1} \left( \log \frac{\rho_i - 1}{\rho_i} \right)^2 \right] \quad (3-9) $$

où dans le deuxième terme $\rho_i=\theta_{ml}$ pour $i=1$, sachant implicitement $\rho_i$ pour $i=1-1$ et $i=1+1$.

[Grandis, 1994] a validé cet algorithme à l'aide de données synthétiques bruitées, et l'a appliqué à l'interprétation de données de terrain [Grandis, 1994; Grandis et al., 1995a; Grandis et al., 1995b]. Les figures 2 et 3 présentent des résultats
<figure><img src="image_112.png" /><figcaption>Figure 2. Exemple d'inversion de données synthétiques correspondant à un milieu tabulaire à l'aide de la méthode proposée par [Grandis, 1994] et [Grandis et al., 1995b]. Les données sont les impédances bruitées à 5%, calculées avec un échantillonnage régulier de 5 points par décade sur une échelle logarithmique de périodes, entre 0,0025 et 250 secondes : (a) résultats sous forme de modèle moyen (trait continu) et modèle déduit de la probabilité a posteriori maximale (pointillé); (b) probabilité a posteriori de la résistivité dans chacune des couches, correspondant à un intervalle de confiance à environ 68%, et modèle synthétique (tireté) (d'après [Grandis, 1994])</figcaption></figure>

obtenus pour des données synthétiques bruitées et des données de terrain.

### 3.2.3. Cas de plaques minces hétérogènes.
L'approximation des plaques minces permet de simplifier la résolution du problème direct de l'induction pour des structures dont l'épaisseur est petite devant leur extension latérale et devant la profondeur de pénétration du champ inducteur. De nombreuses structures rencontrées en géophysique peuvent être décrites comme des plaques minces, et de nombreux auteurs ont utilisé ce formalisme pour décrire le comportement du champ électromagnétique terrestre (cf. par exemple [Weidelt, 1975; Vasseur et Weidelt, 1977; Dawson et Weaver, 1979; Ranganayaki et Madden, 1980; Mareschal et Vasseur, 1984; McKirdy et al., 1985; Robertson, 1987; Tarits, 1989; Terra, 1993; Menvielle et

Tarits, 1994]. Dans le cas présent, la plaque mince est située en surface ou à l'intérieur d'un milieu tabulaire connu; elle comporte une zone hétérogène d'extension limitée. La zone hétérogène est décomposée en N mailles homogènes de conductivité intégrée $\tau = \tau_n + \tau_a$, $\tau_n$ étant la conductivité intégrée de la plaque mince à l'extérieur de la zone hétérogène (Figure 4).

Le champ électrique à la surface de la plaque est solution de :

$$E(M_j) = E_n(M_j) - i\omega\mu_0 \sum_{k=1, N} \sigma_a(M_k) G(M_j, M_k) E(M_k)$$

où $G(M_j, M_k)$ est le tenseur de Green G intégré sur la maille $M_j$. Ce tenseur (2x2) représente les composantes du champ électrique créé en $M_j$ par un dipôle électrique unitaire placé au centre de la
<figure><img src="image_113.png" /><figcaption>Figure 3. Exemple d'inversion de données de terrain. Les résultats sont présentés sous forme de probabilité à posteriori de la résistivité dans chacune des couches, correspondant à un intervalle de confiance à environ 68% (d'après [Grandis, 1994]).</figcaption></figure>

maille $M_k$ de conductance $\tau_n$. La résolution complète du problème direct nécessite donc de résoudre un système linéaire $2Nx2N$ pour calculer le champ électrique à la surface de la plaque mince, puis d'en déduire les observables correspondants en chacun des points de mesure. Ce système est résolu par la méthode itérative de Gauss Seidel, et [Roussignol et al., 1993] ont proposé de l'effectuer qu'un nombre limité d'itérations à chaque étape du balayage, calculant ainsi simultanément le champ électrique et les conductances $\tau_a$ à l'aide d'une chaîne de Markov à deux dimensions.

[Jouanne, 1991] et [Grandis, 1994] ont montré expérimentalement que cette chaîne converge pour les situations rencontrées en géophysique dès l'instant où deux itérations de Gauss-Seidel sont effectuées à chaque étape du balayage. Cela signifie que le calcul approché du champ électrique pour la maille $k_c$ lors de l'actualisation de la maille $k$ ne peut être effectué en négligeant complètement l'effet secondaire produit sur toutes les mailles

différentes de $k$ et de $k_c$ par le changement de conductance sur la maille $k$, effet qui commence à être pris en compte lors de la deuxième itération de Gauss Seidel.

La chaîne de Markov ainsi définie est une chaîne de Markov à deux dimensions. La convergence de cette chaîne a été étudiée expérimentalement par [Grandis, 1994] à l'aide de modèles synthétiques simples, dont un exemple est présenté figure 4. Les observables utilisés ont été le tenseur d'impédance ou le vecteur d'induction, calculé aux différents points d'observation pour différentes périodes et bruité à l'aide d'une variable aléatoire suivant une distribution gaussienne $N(0, \sigma)$. Le résultat de l'inversion est la loi *a posteriori* des paramètres, représentée par les probabilités marginales de chacune des conductances des mailles. La figure 5 présente un exemple des résultats obtenus pour une plaque mince hétérogène située en surface.

Les résultats présentés figure 5 montrent en particulier que l'information contenus dans les vecteurs d'induction permet de
<figure><img src="image_114.png" /></figure>

<table><thead><tr><th></th><th></th><th colspan="2">surface</th></tr></thead><tbody><tr><td>10 : 500 Ohm.m</td><td>5 km</td><td>500 Ohm.m</td><td>10 km</td></tr><tr><td>500 Ohm.m</td><td>40 km</td><td>10 : 500 Ohm.m</td><td>5 km</td></tr><tr><td>1000 Ohm.m</td><td>70 km</td><td>500 Ohm.m</td><td>40 km</td></tr><tr><td>10 Ohm.m</td><td></td><td>1000 Ohm.m</td><td>70 km</td></tr><tr><td></td><td></td><td>10 Ohm.m</td><td></td></tr></tbody></table>

Figure 4. Modèle synthétique utilisé dans le cas de l'inversion de plaques minces hétérogènes. (a) maillage et répartition de conductance; (b) milieu stratifié pour le modèle comportant une plaque mince en surface; (c) milieu normal stratifié pour le modèle comportant une plaque mince en profondeur. La notation '10:500 Ohm.m' représente le contraste de résistivité entre la partie homogène de la plaque mince et l'hétérogénéité résistante. Dans le cas où la plaque mince est résistante avec une hétérogénéité conductrice, le contraste est simplement inversé (500:10 Ohm.m) (d'après [Grandis, 1994]).

déccrire avec une bonne précision les zones anormalement conductrices. Ceci s'explique aisément du point de vue physique. Les vecteurs d'induction caractérisent les variations anormales du champ magnétique engendrées par les courants circulant dans les zones anormalement conductrices du milieu sous-jaçent (cf. par exemple, [Parkinson, 1959; Wiese, 1965]. Lorsque la fréquence du champ primaire inducteur est suffisamment basse, l'intensité de ces courants, et donc celle du champ magnétique associé, est proportionnelle à la conductivité intégrée de la structure anormalement conductrice [Le Mouel et Menvielle, 1982; Menvielle, 1984; Menvielle, 1988] Cette propriété asymptotique a notamment été vérifiée lors de l'étude numérique de l'effet de bord de mer [Tarits, 1982; Menvielle et al., 1982].

Nous avons appliqué cette approche à l'étude de la conductivité intégrée de la croûte en Finlande, ainsi qu'en Tchéquie et en Slovaquie, en utilisant les données de sondages géomagnétiques profond et de sondages magnétotelluriques disponibles dans ces régions. Les résultats obtenus pour la Tchéquie et la Slovaquie, présentés figure 6, mettent en évidence l'intérêt géodynamique de cette cartographie. Le trajet des anomalies de conductivité déduites de l'analyse de la géométrie des vecteurs d'induction est confirmé, et souvent précisé, et la cartographie fait apparaître une réelle corrélation entre la tectonique régionale et la distribution de la conductivité intégrée dans la croûte.

L'inversion confirme en effet la présence d'une anomalie de conductivité sous le Massif des Carpathes. Cette anomalie, déja étudiée par de nombreux auteurs (cf. par exemple, [Praus et al., 1981; Cerv et al., 1984; Jankowski et al., 1985; Cerv et al., 1987], correspond à la subduction des Carpathes extérieures sous les Carpathes intérieures. Mais la carte obtenue met en outre en évidence des anomalies qui coincident avec des accidents tectoniques majeurs situés au Nord des Carpathes, tels que le linéamant Moravo-Silésien et celui de Labe. Ces anomalies, jamais décrites parce que difficiles à mettre en évidence par

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

Figure 5. Résultat de l'inversion dans le cas d'une plaque mince en profondeur avec une répartition logarithmique des conductances a priori θm (m=1, ... M; voir texte section 3.2.1 pour plus de détails) Le modèle est décrit Figure 4. Les données inversées sont calculés au centre de chaque maille de la zone décrite Figure 4a et bruitées à 10%. Les probabilités marginales a posteriori de conductance dans chaque maille sont représentées par la probabilité de transition à la vingtième itération. (a) inversion du tenseur d'impédance, hétérogénéité résistante; (b) inversion du tenseur d'impédance, hétérogénéité conductrice; (c) inversion du vecteur d'induction, hétérogénéité résistante; (d) inversion du vecteur d'induction, hétérogénéité conductrice.
<figure><img src="image_116.png" /><figcaption>Figure 6.a. Exemple de cartographie de la conductivité électrique de la lithosphère à l'aide de l'algorithme basé sur l'approximation des plaques minces. Les données inversées sont les vecteurs d'induction observés pour une période de 20 minutes (a : partie réelle; b : partie imaginaire);</figcaption></figure>

l'interprétation qualitative des cartes de vecteur d'induction, constituent un apport original des méthodes d'imagerie électromagnétique.

Ces résultats montrent que le traitement de cartes de vecteurs d'induction à l'aide de méthodes quantitatives d'imagerie électromagnétique fournit des informations originales sur la structure de la croûte. A la suite des nombreuses campagnes de sondages géomagnétiques effectuées depuis le début des années soixante, des cartes de vecteur d'induction ont été publiées pour de nombreuses régions du monde. La réinterprétation quantitative de ces vecteurs d'induction sous forme de cartes de la conductivité intégrée de la croûte à l'échelle régionale constituera très vraisemblablement une contribution majeure de l'électromagnétisme à la tectonique et la géodynamique.

### 3.2.4. Cas bi-dimensionnel.
L'approximation bi-dimensionnelle (2-D) permet de simplifier la résolution du problème direct de l'induction pour des structures dont l'extension dans une des directions du plan horizontal est grande devant leur extension latérale et devant la profondeur de pénétration du champ inducteur. La structure est alors définie par

la distribution de résistivité dans un plan $Σ$ perpendiculaire à son axe $Δ$. En situation 2-D, le champ électromagnétique se découle en deux modes séparés, correspondant l'un à un champ électrique horizontal parallèle à $Δ$ et à un champ magnétique perpendiculaire à $Δ$ (mode $E_{||}$), et l'autre à un champ électrique horizontal parallèle à $Δ$ et à un champ électrique perpendiculaire à $Δ$ (mode $B_{||}$). Ce découplage simplifie relativement les calculs, et de nombreux auteurs ont utilisé ce formalisme pour décrire le comportement du champ électromagnétique terrestre au voisinage de structures tectoniques allongées.

La résolution du problème direct pour des structures 2-D est maintenant classique, mais elle reste cependant très coûteuse en temps calcul. Les méthodes inverses qui ont été développées pour des structures 2-D utilisent donc en général une solution approchée du problème direct, dans le cadre d'hypothèses sur la loi *a priori* des paramètres (cf. par exemple [Mackie et al., 1988] ou sur le caractère lisse du modèle (cf. par exemple [Smith et Booker, 1988]). Dans cette perspective, la solution la plus simple est celle qui consiste à juxtaposer les profils 1-D déduits des observations faites en chacune des stations, et à prendre pour
<figure><img src="image_117.png" /><figcaption>Figure 6.b Modèle moyen a posteriori après lissage sur une réseau de mailles deux tois plus dense que celui utilisé lors de l'inversion.LL:Linéament de Labe; MSL : linéament Moravo-Silésien; BMA : anomalie de conductivité du Massif Bohémien; CPA : anomalie de conductivité des Carpathes (d'après [Grandis, 1994]).</figcaption></figure>

solution du problème inverse le résultat du lissage de la distri-
bution de résistivité ainsi obtenue.

La traduction de cette démarche sous forme d'un algorithme
d'inversion bidimensionnel est aisée dans le cadre de l'approche
développée par [Roussignol et al., 1993] et [Grandis, 1994] :

-le domaine de Σ dans lequel on recherche la distribution de résistivité est décrit à l'aide de mailles rectangulaires homogènes : les paramètres du modèle sont donc les résistivités de chacune des mailles;

-l'algorithme d'inversion est analoque à celui utilisé dans le cas 1-D, présenté section 3.2.2, mais avec une loi a priori pénalisant des sauts de résistivité entre deux mailles voisines dans Σ;

-les données à inverser sont les impédances associées au mode E_{||}. Le champ électromagnétique correspondant à ce mode est en effet caractérisé par une distance d'ajustement latérale faible, c'est à dire qu'il n'est perturbé par un contraste de conductivité que dans le voisinage immédiat de ce contraste. Le champ électromagnétique correspondant au mode B_{||} est, au contraire, caractérisé par une distance d'ajustement latérale pouvant atteindre plusieurs centaines, voire plusieurs milliers, de kilomètres [Ranganayaki et Madden, 1980; Dawson et al., 1982; Menvielle, 1988; Weaver et Dawson, 1992; Menvielle et Tarits, 1994]. En ce qui concerne les impédances, cela signifie que l'impédance associée au mode E_{||} n'est pas sensible à la présence d'éventuelles structures situées loin du domaine étudié, alors que celle associée au mode B_{||} peut y être sensible.

[Touijar, 1994] et [Menvielle et Roussignol, 1995] ont développé un tel algorithme d'inversion pour des structures bidimensionnelles, et l'ont validé pour des données synthétiques. La figure 7 présente un exemple des résultats obtenus.

4. CONCLUSION

Ainsi que le montrent les résultats présentés dans cet article, la
contribution des équipes françaises constituent une étape dans le
développement de méthodes d'imagerie électromagnétique. Les
méthodes développées par ces équipes permettent en effet, d'une
part, de caractériser l'ensemble des modèles acceptables et,
d'autre part, de quantifier l'information apportée par les données.
Les résultats obtenus pour des structures de géométrie variée
justifient les choix effectués, et notamment l'utilisation de
chaînes de Markov ergodiques dont la loi de transition se
stabilise asymptotiquement vers la loi a posteriori des
paramètres. En particulier, [Grandis, 1994] a montré que la
chaîne de Markov proposée par [Roussignol et al., 1993]
converge rapidement sur des exemples représentatifs de
situations géophysiques, ce qui permet d'estimer la loi a
posteriori des paramètres à des coûts en temps calcul qui soient
réalistes. Les méthodes actuellement développées sont cependant
adaptées à l'étude de structures de géométrie particulière : milieu
tabulaire, plaque mince hétérogène dans un milieu tabulaire,
structures bi-dimensionnelles. Tel n'est pas le cas des structures
rencontrées dans de nombreuses situations géophysiques, et
[Spichak et al., 1995] ont développé une méthode d'imagerie
électromagnétique basée sur l'approche adoptée par [Roussignol
et al., 1993] et permettant l'étude de structures de géométrie
quelconque.

D'un point de vue plus général, les résultats présentés dans cet
article illustrent l'intérêt de l'approche Bayesienne dans la
résolution de problèmes inverses non linéaires en géophysique :
le savoir avant inversion et les hypothèses faites pour
sélectionner les modèles retenus sont objectivés et décrits quanti-
tativement par la loi a priori ; la résolution du problème inverse
a

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

b

<figure><img src="image_119.png" /><figcaption>Figure 7. Exemple d'inversion d'un modèle synthétique bidimensionnel. (a) Le modèle; (b) le résultat de l'inversion des impédances bruitées correspondant au mode $E_{||}$, calculées pour 26 périodes comprises entre 0,0025 et 250 secondes, en 9 points équirépartis le long d'un profil au dessus du domaine représenté. La largeur de ce domaine est de 90 km, et sa hauteur de 20 km.</figcaption></figure>
est possible pour une grande variété de lois *a priori*, sans qu'il soit nécessaire de linéariser le problème direct; l'information disponible sur les paramètres du modèle après inversion est décrite quantitativement par la loi *a posteriori*.

REFERENCES

Andersen, R.S. et Seneta, E. A simple statistical estimation procedure for Monte Carlo inversion in geophysics, *PAGEOPH*, *91*, 5-13, 1971.

Backus, G. Bayesian inference in geomagnetism, *Geophys. J.*, **92**, 125-142, 1988.

Bayes, T. An essay toward solving a problem in the doctrine of chances, *Phil. Trans. R. Soc. Lond.*, **53**, 370-418, 1763.

Berger, J.O. Statistical decision theory and Bayesian analysis, *Springer-Verlag, New York*, 1985.

Box, G., Tiao, G. Bayesian inference in statistical analysis, *Addison-Wesley, Reading*, 1973.

Cerv, V., Pek, J. et Praus, O. Models of geoelectrical anomalies in Czechoslovakia, *J. Geophys.*, **55**, 161-168, 1984.

Cerv, V., Pek, J. et Praus, O. Numerical modelling of geoelectrical structures in Czechoslovakia, *Phys. Earth Planet. Inter.*, **45**, 170-178, 1987.

Constable, S.C., Parker, R.L. et Constable, C.G. Occam's inversion, a practical algorithm for generating smooth models from electromagnetic sounding data, *Geophysics*, **52**, 289-300, 1987.

Dawson, T.W., Weaver, J.T. et Raval, U. B-polarization induction in two generalized thin sheets at the surface of a conducting half space, *Geophys. J. Roy. astr. Soc.*, **69**, 209-234, 1982.

Dawson, T.W. et Weaver, J.T. Three dimensional induction in non uniform thin sheet at the surface of a uniformly conducting earth, *Geophys. J. Roy. astr. Soc.*, **59**, 445-462, 1979.

Dosso, S.E. et Oldenburg, D.W. Magnetotelluric appraisal using simulated annealing, *Geophys. J. Int.*, **106**, 379-385, 1991.

Everett, M.E. et Schultz, A. Two-dimensional non linear magnetotelluric inversion using a genetic algorithm, *J. Geomag. Geoelectr*, **45**, 1013-1026, 1993.

Filloux, J.H. Ocean floor magnetotelluric sounding over North Central Pacific, *Nature*, **269**, 297-301, 1977.

Filloux, J.H. Magnetotelluric soundings over the north-west Pacific may reveal spatial dependence of depth and conductance of the asthenosphere, *Earth Planet. Sci. Lett.*, **46**, 244-252, 1980.

Fisher, G. et Le Quang, B.V. Topography and minimization of the standard deviation in one dimensional magnetotelluric modelling, *Geophys. J. Roy. astr. Soc.*, **67**, 279-292, 1981.

Gibert, D., Tournerie, B. et Virieux, J. High resolution electromagnetic imaging of the conductive Earth interior, **10**, 341-351, 1994.

Gibert, D. et Virieux, J. Electromagnetic imaging and simulated annealing, *J. Geophys. Res.*, **96**, 8057-8067, 1991.

Grandis, H. Imagerie electromagnetique Bayesienne par la simulation d'une chaîne de Markov, *Doctorat d'Université. Université Paris VII*, 1994.

Grandis, H., Jouanne, V., Dupis, A., Menvielle, M. et Tarits, P. Etude électromagnétique de la partie Nord du horst de Couy (Cher, France). *Bull. Soc. géol. France, soumis*, 1995a.

Grandis, H., Menvielle, M. et Roussignol, M. Bayesian electromagnetic imaging with Markov chains: example of the 1-D magnetotelluric inverse problem. *Geophys. J. Int.*, *soumis*, 1995b.

Jankowski, J., Tarlowski, Z., Praus, O., Pecova, J. et Petr, V. The results of deep geomagnetic soundings in the west Carpathian, *Geophys. J. Roy. astr. Soc.*, **80**, 561-574, 1985.

Jouanne, V., Menvielle, M., Roussignol, M. et Tarits, P. 3-D Bayesian electromagnetic imaging, communication orale, Assemblée Générale de l'E.G.S. Barcelone, 1989.

Jouanne, V. Application des techniques statistiques Bayesiennes à l'inversion de données électromagnétiques, *Doctorat d'Université. Université Paris VII*, 1991.

Law, L.K. et Greenhouse, J.P. Geomagnetic variation soundings of the asthenosphere beneath the Juan de Fuca Ridge, *J. Geophys. Res.*, **86**, 967-978, 1981.

Le Mouel, J.L. et Menvielle, M. Geomagnetic variation anomalies and deflection of telluric currents, *Geophys. J. Roy. astr. Soc.*, **68**, 575-587, 1982.

Mackie, R.L., Bennett, B.R. et Madden, T.R. Long period magnetotelluric measurements near the central California coast: a land-locked view of the conductivity structure under the Pacific Ocean, *Geophys. J.*, **95**, 181-194, 1988.

Mareschal, M. et Vasseur, G. Bimodal induction in non-uniform thin sheets: do the presents algorithms work for regional studies? *J. Geophys.*, **55**, 203-213, 1984.

Marroquin, J., Mitter, S. et Poggio, T. Probabilistic solution of ill-posed problems in computational vision, *J. Amer. Statist. Assoc.*, **82**, 76-89, 1987.

McKirdy, D.M., Weaver, J.T. et Dawson, T.W. Induction in a thin sheet of variable conductance at the surface of a stratified earth -II: three dimensional theory, *Geophys. J. Roy. astr. Soc.*, **80**, 177-194, 1985.

Menke, W. Geophysical data analysis: discrete inverse theory, *Academic Press*, 1984.

Menvielle, M., Rossignol, J.C. et Tarits, P. The coast effect in terms of deviated electric currents: a numerical study, *Phys. Earth Planet. Inter.*, **28**, 118-128, 1982.

Menvielle, M. Anomalies de conductivité et phénomènes de canalisation, *Doctorat d'Etat. Université Paris VI*, 1984.

Menvielle, M. Effects of crustal conductivity heterogeneities on the electromagnetic field, *Surv. Geophys.*, **9**, 319-348, 1988.

Menvielle, M. et Roussignol, M. A simple 2-D Bayesian electromagnetic inversion scheme. *Geophys. J. Int.*, *soumis*, 1995.

Menvielle, M. et Tarits, P. Adjustment of the electromagnetic field distorted by 3-D heterogeneities, *Geophys. J. Int.*, **116**, 562-570, 1994.

Metropolis, N., Rosenbluth, A.W., Rosenbluth, M.N., Teller, A.H. et Teller, E. Equation of state calculation by fast computing machines, *J. Chem. Phys.*, **1**, 1087-1092, 1953.

Metropolis, N. et Ulam, S.M. The Monte Carlo method, *J. Amer. Statist. Assoc.*, **44**, 335-341, 1949.

Mosegaard, K. et Tarantola, A. Monte Carlo sampling of solutions to inverse problems. *J. Geophys. Res.*, *soumis*, 1995.

Olhoeff, G.R. Electrical properties of rocks, *Thoulaikian, Y.S., Judd, W.R. et Roy, R.F. editeurs, McGraw Hill*, 257-330, 1980.

Parker, R.L. The inverse problem of electromagnetic induction: existence and construction of solutions based on incomplete data, *J. Geophys. Res.*, **85**, 4421-4428, 1980.

Parkhomenko, E.I. Electrical resistivity of minerals and rocks at high temperature and pressure, *Rev. Geophys. Space Phys.*, **29**, 193-218, 1982.

Parkinson, W.D. Direction of rapid geomagnetic fluctuations, *Geophys. J. Roy. astr. Soc.*, **2**, 1-13, 1959.

Pedersen, J. et Hermance, J.F. Least squares inversion of one dimensional magnetotelluric data: an assessment of procedures employed by Brown University, *Surv. Geophys.*, **8**, 187-231, 1986.

Pedersen, L.B. et Rasmussen, T.M. Inversion of magnetotelluric data: a non linear least squares approach, *Geophys. Prosp.*, **37**, 669-695, 1989.

Praus, O., Pecova, J., Petr, V., Pec, K., Hvozdara, M., Cerv, V., Pek, J. et Lastovickova, L. Electromagnetic induction and electrical conductivity in the earth's body, *Geophysical synthesis in Czechoslovakia, Zatopek, A. editeur, Veda, Bratislava*, 1981.

Press, S.J. Bayesian statistics: principles, models, and application, *John Wiley and Sons*, 1989.

Press, S.J., Flannery, B.P., Teukolsky, S.A. et Vetterling, W.T. Numerical recipes: the art of scientific computing, *Cambridge University Press*, 1989.

Ranganayaki, R. et Madden, T.R. Generalized thin sheet analysis in magnetotellurics: an extension of the Price's analysis, *Geophys. J. Roy. astr. Soc.*, **60**, 445-457, 1980.

Robertson, R.C. The electromagnetic response of thin sheets buried in a uniformly conducting half space, *Geophysics*, **52**, 108-117, 1987.

Roussignol, M., Jouanne, V., Menvielle, M. et Tarits, P. Bayesian electromagnetic imaging, Computer intensive methods, *Hardle, W., Siman, L. editeurs, Physical Verlag*, 85-97, 1993.

Sambridge, M. et Drijonningen, G. Genetic algorithms in seismic waveform inversion, *Geophys. J. Int.*, **109**, 323-342, 1992.

Smith, J.T. et Booker, J.R. Magnetotelluric inversion for minimum structures, *Geophysics*, **53**, 1565-1576, 1988.

Spichak, V. EM-field transformations and their use in interpretation, *Surv. Geophys.*, **11**, 271-301, 1990.

Spichak, V., Menvielle, M. et Roussignol, M. 3-D electromagnetic bayesian imaging. *Proceedings of the G. Hohmann symposium, soumis*, 1995.

Tarantola, A. Inverse problem theory: method for data fitting and model parameters estimation, *Elsevier*, 1987.

Tarits, P. Etude du champ électromagnétique terrestre d'origine intralithosphérique, *Doctorat de Troisième Cycle. Université Paris VII*, 1982.

Tarits, P. Contribution des sondages électromagnétiques profonds à l'étude du manteau supérieur terrestre, *Doctorat d'Etat. Université Paris VII*, 1989.

Tarits, P., Jouanne, V., Menvielle, M. et Roussignol, M. Bayesian statistics of non-linear inverse problems: example of the magnetotelluric 1-D inverse problem, *Geophys. J. Int.*, **119**, 353-368, 1994.

Tarits, P. et Jouanne, V. Résultats de sondages magnétotelluriques sous-marins et structure thermique des points chauds, *Bull. Soc. géol.*
France, VI, 921-931, 1990.

**Terra, A.** Modélisation électromagnétique de deux plaques minces hétérogènes: application à l'étude de la croûte et du manteau, *Doctorat d'Université. Université Paris VII*, 1993.

**Tikhonov, A.N., Arsenin, V.Y.** Solution of ill-posed problems, *John Wiley and Sons*, 1977.

**Touijar, D.** Algorithmes stochastiques et application à l'imagerie électromagnétique Bayesienne, *Doctorat d'Université. Université de Lille*, 1994.

**Tournerie, B., Gibert, D. et Virieux, J.** Inversion of the COPROD2 magnetotelluric data using a diffusive to propagative mapping. *Geophys. Res. Lett., soumis*, 1995.

**Vasseur, G. et Weidelt, P.** Bimodal electromagnetic induction in non-uniform thin sheets with an application to the Northern Pyrenean induction anomaly, *Geophys. J. Roy. astr. Soc.*, **51**, 669-690, 1977.

**Virieux, J., Flores-Luna, C. et Gibert, D.** Asymptotic theory for diffusive electromagnetic imaging, *Geophys. J. Int.*, **119**, 857-868, 1994.

**Weaver, J.T. et Dawson, T.W.** Adjustment distance in TM mode electromagnetic induction, *Geophys. J. Int.*, **108**, 293-300, 1992.

**Weidelt, P.** Electromagnetic induction in three-dimensional structures, *J. Geophys.*, **41**, 85-109, 1975.

**Wiese, H.** Geomagnetische induktionspfeile in der CSSR hervorgeraufendurch grossräumige elektrische lietfähigkeitsstrukturen, *Studia geophys. geod.*, **9**, 415-419, 1965.

**Wiggins, R.A.** Monte Carlo inversion of body wave observations, *J. Geophys. Res.*, **74**, 3171-3181, 1969.

**Wiggins, R.A.** The general linear inverse problem: implication of surface waves and free oscillations for earth structure, *Rev. Geophys. Space Phys.*, **10**, 251-285, 1972.

**Wu, F.T.** The inverse problem of magnetotelluric sounding, *Geophysics*, **33**, 972-979, 1968.

**Yukutake, T., Filloux, J.H., Segawa, J., Hamano, Y. et Utada, H.** Preliminary report on a magnetotelluric array study in the northwest Pacific, *J. Geomag. Geoelectr.*, **35**, 575-587, 1983.
# LE CHAMP GEOMAGNETIQUE ET LA CONVECTION DANS LE NOYAU

Gauthier HULOT et Dominique JAULT

Laboratoire de Géomagnétisme, Institut de Physique du Globe de Paris, 4 Place Jussieu, Tour 24-25, 75252 Paris Cedex 05

Nous présentons des résultats récents concernant le champ géomagnétique principal, en insistant plus particulièrement sur la contribution française (ENS de Paris, IPGP, EOPGS) durant les quatre dernières années. Une première série de travaux portant sur l'analyse du champ géomagnétique observé fournit d'abord l'occasion de rappeler les caractéristiques principales de ce champ. Quelques considérations sur la dynamique du noyau liquide de la Terre permettent alors de souligner ce qui fait la spécificité de la géodynamo. Un certain nombre de résultats récemment obtenus par modélisation numérique sont ensuite présentés et permettent d'apprécier les conditions physiques dans lesquels la géodynamo fonctionne vraisemblablement. On peut alors interpréter les observations, variation séculaire et secousses géomagnétiques, en termes de mouvements animant le fluide à la surface du noyau, de telle manière qu'il est même possible d'en étudier les liens avec les variations décennales de la rotation de la Terre (irrégularités de la durée du jour et mouvements du pôle de rotation).

We present recent results concerning the Earth's main magnetic field with special emphasis on the french contribution (ENS de Paris, IPGP, EOPGS) during the last four years. The presentation of a series of studies dealing with the analysis of the observed geomagnetic field first provides an opportunity to recall the main characteristics of this field. Some considerations on the dynamics of the Earth's liquid core then makes it possible to underline specific features of the geodynamo. Results recently derived from numerical modelling are also presented and allow some inference to be made about the physical conditions under which the geodynamo is likely to operate. Eventually geomagnetic observations (i.e., secular variation and geomagnetic jerks) are interpreted in terms of flows acting at the top of the liquid core, in such a way that some connections can also be made with the decade variations of the Earth's rotation (length of day variations and motion of the pole of rotation).

## 1. GRANDS TRAITS DU CHAMP GEOMAGNETIQUE CONTEMPORAIN

### 1.1. Modélisation du champ magnétique principal et de sa variation séculaire

Le champ magnétique que nous observons aujourd'hui à la surface de la Terre semble simple en première approximation. Nous pouvons le décrire comme étant principalement d'origine interne, et très semblable au champ qu'engendrerait un dipôle magnétique situé au centre de la Terre, d'intensité environ 8.10²² Am², et dont le pôle Nord (intersection du dipôle avec la surface de la Terre) se situerait au voisinage du point de latitude 79° Sud et de longitude 109° Est (le pôle Nord magnétique est donc dans l'hémisphère Sud). Cette apparente simplicité cache cependant une réalité plus complexe. Ainsi le champ magnétique possède également des sources dans l'ionosphère et la magnétosphère. Ces sources, dites externes, sont de nature très différente des sources internes et ne seront pas considérées au cours de cet article. Leur modeste contribution, le champ magnétique d'origine externe, peut être caractérisée et soustraite des observations [e.g., Gire, 1985; Cohen, 1989], ce qui permet d'accéder au champ magnétique d'origine interne que l'on décrit traditionnellement au moyen d'un développement en harmonique sphérique:

$$ \vec{B} = -\vec{V}V; \\ V = a \sum_{n=1}^{\infty} \sum_{m=0}^{n} (g_n^m Y_n^{mc}(\theta, \phi) + h_n^m Y_n^{ms}(\theta, \phi)) \left(\frac{a}{r}\right)^{n+1} \quad (1) $$

où a est le rayon de la Terre, r, θ, φ sont les coordonnées sphériques et les $Y_n^m(c, s)$ sont les harmoniques élémentaires de degré n et ordre m. Le degré 1 décrit la composante dipolaire du champ, et les degrés plus élevés correspondent à des structures d'échelle spatiale d'autant plus petite que n et m sont grands.

La donnée d'un jeu de coefficients de Gauss $g_n^m$ et $h_n^m$ définit complètement le champ d'origine interne et constitue ce que l'on

appelle un "modèle géomagnétique". L'étude de modèles très détaillés obtenus grâce au satellite MAGSAT pour l'époque 1980 permet de mettre en évidence l'existence, au sein du champ d'origine interne, de deux contributions d'origines physiques et géographiques très distinctes: un champ magnétique "crustal", dont les sources sont des roches aimantées situées dans les quarante à cinquante premiers kilomètres de la Terre solide (profondeur de l'isotherme au delà de laquelle tout ferromagnétisme devient impossible et qui correspond, grosso modo, à l'épaisseur de la croûte terrestre), et un champ magnétique "nucléaire" dont les sources sont des courants électriques parcourant le noyau liquide et conducteur de la Terre. Ces deux champs possèdent des comportements spatiaux très différents à cause principalement de la distance relative des sources crustales, proches de la surface de la Terre, et des sources nucléaires, plus lointaines. Ceci se traduit par le fait que chaque degré n (soit en gros chaque échelle spatiale) du champ crustal contribue à peu près également à la valeur quadratique moyenne

Le champ magnétique principal ne se caractérise pas seulement par une structure spatiale complexe. Il évolue sur des échelles de temps allant typiquement de quelques années à, pour ainsi dire, l'âge de la Terre [e.g., Courtillot et Le Mouël, 1988; Prévot et Perrin, 1992]. Malheureusement, et pour d'évidentes raisons
<figure><img src="image_120.png" /></figure>

Figure 1. Cartes de la composante radiale du champ magnétique principal à la base du manteau en 1700 (a), 1840 (b), 1980 (c). (D'après Bloxham et Jackson, 1992). Isocontours tous les 100 µT. Les repères 1 et 2 situent deux exemples de structures apparemment stationnaires. Par ailleurs, on a représenté en grisé les taches définies par des isocontours $B_r = 0$. Celle associée à une variation apparemment forte du flux se trouve à l'aplomb de l'Afrique du Sud. Cette variation est récente et est responsable de l'accroissement brutal (par fusion avec la "langue" descendant de l'équateur le long de l'Amérique du Sud), de la grande tache positive dans l'hémisphère Sud.

pratiques, il n'est pas possible de suivre cette évolution à travers les temps de manière homogène. On distingue ainsi la "variation séculaire", directement observable à l'échelle historique, des variations archéomagnétiques et paléomagnétiques, indirec-tement observables par le biais de matériaux archéologiques et géologiques aimantés dans le passé. Dans la mesure où seule l'observation directe permet, pour l'instant, d'obtenir suffisamment de données précisément datées pour autoriser la construction d'un modèle géomagnétique variant dans le temps, c'est avant tout à la variation séculaire, aux informations qu'elle contient, et au mécanisme qui en est à l'origine, que nous allons nous intéresser. Il importera cependant de replacer les conclusions auxquelles nous arriverons dans le contexte plus vaste d'un champ géomagnétique vieux de plus de trois milliard d'années.

Si la couverture géographique exceptionnelle des mesures faites par le satellite MAGSAT a permis de déterminer le champ magnétique principal avec précision, la relative brièveté de son vol n'a pas permis d'améliorer la détermination de la variation séculaire. Celle-ci repose depuis près de trois siècles sur des mesures réalisées au cours de campagnes (sur Terre ou en mer) ainsi que sur des enregistrements continues effectués dans des observatoires (au nombre de 150 environ aujourd'hui). Complétant le modèle de champ magnétique principal pour l'époque 1980 {$g_n^m, h_n^m$} (limité au degré 13 du fait de la présence du champ crustal), nous disposons donc de modèles de variation séculaire {$g_n^{\tilde{m}}, h_n^{\tilde{m}}$}, moins fiables certes (jusqu'au degré 7 environ), mais variant dans le temps. En pratique et en général, nous utiliserons un modèle géomagnétique variant sur l'ensemble des trois siècles, {$g_n^m(t), h_n^m(t)$}, publié par Bloxham et Jackson [1992].

## 1.2. Quelle stationnarité pour le champ magnétique principal ?

Disposant ainsi d'un "film" du champ magnétique principal et utilisant un prolongement à travers le manteau considéré comme isolant (cf. Hulot [1992] et Poirier et Le Mouël [1992], pour des discussions récentes de la validité de cette hypothèse), il est alors possible de l'étudier là où il est généré, à la base du manteau (Figure 1). Plusieurs types d'études peuvent être, et ont été, entreprises. Bloxham et al. [1989] ont ainsi soigneusement analysé l'évolution de structures caractéristiques (telles que les courbes $\Gamma_{B_r=0}$ d'isovaleurs $B_r = 0$), dans le double but de retrouver des propriétés attendues d'un point de vue théorique (comme la conservation du flux magnétique à travers les surfaces délimitées par ces courbes $\Gamma_{B_r=0}$, nous y reviendrons ultérieurement) et de distinguer les structures stationnaires des structures variables. On constate en effet que certaines de ces structures semblent se maintenir pendant toute la durée d'observation tandis que le reste du champ évolue substantiellement (Figure 1). Se pose alors la question de la portée d'une telle observation. Combien de temps peut-on ainsi imaginer que de telles taches puissent se maintenir ? Certains auteurs proposent de les associer au champ dipolaire axial (que l'on sait de polarité stable sur des dizaines de milliers d'années) et pensent que l'on peut en retrouver la trace, très atténuée cependant, jusque dans les enregistrements paléomagnétiques des quelques derniers Ma [e.g., Gubbins et Kelly, 1994]. Mais il importe d'abord d'évaluer le caractère remarquable de la persistance de ces taches. Ces dernières ne constituent-elles pas tout simplement une manifestation temporaire normale du champ non dipôle, et ne sont elles pas amenées à se métamorphoser dans les quelques siècles à venir ?

Des éléments de réponse peuvent être obtenus en adoptant un point de vue statistique. Ainsi Hulot et Le Mouël [1994] ont récemment montré que la variabilité du champ non dipôle est statistiquement semblable à celle qu'engendrerait un processus stochastique relativement simple. Chaque coefficient de Gauss (par exemple $g_n^m(t)$) peut ainsi être indépendamment considéré comme une fonction aléatoire caractérisée par une moyenne nulle, une variance $\sigma_n$, et un temps de corrélation typique $\tau_c(n)$. Les spectres du champ magnétique principal et de sa dérivée temporelle sont alors susceptibles de fluctuer autour de spectres moyens $W_\sigma(n)$ et $W_\sigma^*(n)$, lesquels suivent chacun une loi exponentielle en fonction du degré n (Figure 2). L'ensemble du processus stochastique est alors caractérisé par ces deux spectres qui, notons-le bien, peuvent être déterminés
<figure><img src="image_121.png" /></figure>

Figure 2. Spectres observés en 1980 pour la Variation séculaire (a, en (nTan⁻¹)²), et pour le champ principal (b, (nT)²) à la base du manteau. Les carré noirs représentent les intervalles à l'intérieur desquelles doivent passer les spectres théoriques, compte tenu des fluctuations statistiques attendues. Les différentes droites représentent des spectres compatibles avec les données. On notera l'échelle logarithmique. D'après Hulot et Le Mouël [1994].

dès que l'on dispose de quelques années d'observations. En particulier la relation $\tau_c(n) = (W_\sigma(n)/W_\sigma(n))^{1/2}$ définit complètement le temps de corrélation typique du champ non dipôle de degré n. Ce temps $\tau_c(n)$ décroît d'environ 200 à environ 100 ans lorsque l'on passe du degré 2 au degré 7. L'étude de Hulot et Le Mouël [1994] permet alors de conclure que les taches vues par Bloxham et al. [1989] n'ont pas encore été observées assez longtemps pour pouvoir être considérées comme étant réellement stationnaires.

En fait il apparaît que la durée minimum d'observation permettant de séparer les parties stationnaires et variables du champ principal devrait être de l'ordre du millier d'années, voir plus si l'on s'intéresse aussi à la variation non-axiale du champ dipôle. C'est pourquoi il nous faut maintenant nous tourner vers l'archéomagnétisme et compléter l'étude précédente. Nous disposons déjà d'une bonne base de données récemment réactualisée pour la France [Bucur, 1994] et pour le reste du monde [Daly et Le Goff, 1995]. Une analyse rapide de cette base de données suggère déjà que la composante équatoriale du champ dipôle se conformerait lui aussi au modèle stochastique. Une analyse plus fine reste cependant à mener.

Ainsi la variation "séculaire", dont le nom même rappelle qu'il s'agit de la variation que l'on observe en pratique depuis quelques siècles, possède précisément des constantes de temps séculaires et semble représenter l'essentiel des variations du champ magnétique principal, à l'exclusion du champ dipôle axial qui est dominé par des variations à beaucoup plus long terme (constantes de temps des inversions). Cette vision simple d'un champ fluctuant rapidement autour d'un champ moyen beaucoup

plus lentement variable, semble d'ailleurs bien s'accorder avec les observations paléomagnétiques [Courtillot et al., 1992], même si l'on sait que le champ moyen ne peut être rigoureusement dipolaire axial [e.g., Quidelleur et al., 1994; Gubbins et Kelly, 1994]. Elle permet en tout cas déjà de jeter des ponts entre géomagnétisme et paléomagnétisme [e.g., Hulot et al., 1994].

### 1.3. Secousses géomagnétiques

Si l'étude précédente met l'accent sur le fait que l'évolution du champ magnétique principal est en première approximation régulière à l'échelle séculaire, elle ne doit pas cacher le fait que l'on observe parfois des "secousses" géomagnétiques à l'échelle planétaire. Ces secousses ("jerks" en anglais) interviennent à intervalles irréguliers et se manifestent sous la forme d'un changement de tendance, très rapide (en quelques années) et quasi-simultanée, de la variation séculaire en tout point du Globe (Figure 3). Mises en évidence depuis relativement peu de temps [Courtillot et al., 1978], ces secousses ont maintenant été révélées dans un grand nombre d'enregistrement d'observatoires [e.g., Courtillot et Le Mouël, 1988]. Si diverses analyses en termes d'harmoniques sphériques ont permis de les caractériser en partie [e.g., Stewart, 1991; McLeod, 1992], c'est l'étude toute récente de Alexandrescu et al. [1995] qui a permis de lever les derniers doutes [Alldredge, 1984; Backus et al., 1987] quant à la réalité et l'origine interne du phénomène (Figure 3). Utilisant une analyse en ondelettes originale, cette étude révèle directement les secousses et donne en outre accès à des informations nouvelles sur leur nature (il ne s'agit pas, comme on le croyait jusqu'alors, d'un simple "saut" de la dérivée temporelle seconde du champ magnétique principal). Ces informations devraient nous permettre de préciser les mécanismes qui sont à l'origine des secousses, et de renouveler les contraintes qui portent sur la faiblesse de la conductivité du manteau [e.g., Achache et al., 1980; Backus, 1983].

## 2. DYNAMIQUE DU NOYAU LIQUIDE

L'existence et l'évolution du champ magnétique terrestre sont intimement liées aux mouvements animant le noyau liquide. Au fur et à mesure que la Terre se refroidit, le rayon de la graine solide s'accroît et du fer lourd s'accrète à sa surface, libérant de l'énergie potentielle de gravité et de la chaleur latente de solidification. La partie liquide du noyau se trouve alors brassée par des mouvements de convection qui peuvent entretenir une dynamo, parce que le fluide est bon conducteur d'électricité. Différents groupes de chercheurs, dont un français (collaboration IPGP-ENS Paris), s'efforcent de préciser ce mécanisme.

Parce que la diffusivité magnétique (nombre construit à partir de la conductivité électrique) est très grande devant les autres diffusivités, que ce soit la viscosité cinématique ou les diffusivités chimiques et thermiques, l'étude de la dynamique du noyau doit prendre en compte une gamme d'échelles de temps et d'espace très large. Les profils radiaux de température (avec un écart entre le haut et le bas du noyau liquide qui se compte en milliers de degrés) et de concentration d'éléments légers se mettent en place sur des temps très longs, de l'ordre de l'âge de la Terre, tandis que les mouvements responsables de l'entretien de la dynamo évoluent sur un temps caractéristique beaucoup plus court, de l'ordre de $10^2$ à $10^4$ ans, temps se formant à partir de la diffusivité magnétique. Les études sur l'évolution lente de la température débutent tout juste, mais Le Mouël, Poirier, et Assoumani [1994] ont déjà illustré le rôle central joué par le couplage thermique avec le manteau, moins bon conducteur de la chaleur. Ceci dit, nous ne nous intéresserons ici qu'aux
<figure><img src="image_122.png" /></figure>

Figure 3: Variations du champ magnétique observés sur la composante Y à Chambon-la-Forêt (a), lignes de maxima de la carte temps-fréquence de la transformée en ondelettes du signal précédent (b), variation de la valeur absolue de la transformée en ondelettes le long des lignes de maxima précédentes (c). L'identification d'un jerk se fait grâce à la présence d'un segment rectiligne depuis les hautes fréquences jusqu'aux basses fréquences. C'est ce qu'on observe en (c) pour les lignes de maxima 1, 2, 4, 13 et 14 de (b) correspondant à des jerks en 1901, 1925, 1969 et 1978 pour le signal (a). D'après Alexandrescu et al. [1995].

constantes de temps pertinentes pour le mécanisme de la dynamo proprement dite. L'estimation ci-dessus de ce temps caractéristique montre d'ailleurs que la variation séculaire décrite précédemment est bien la signature des mouvements qui permettent au dipôle de terrestre de se maintenir face à l'effet Joule. Les anomalies de densité et de température évoluant sur cette échelle de temps rapide sont infimes (respectivement de l'ordre de $10^{-5}$ g/cm$^3$ et du millième de degré). Enfin, le champ magnétique est bien régulier parce que toute variation spatiale

brusque se trouverait éliminée par la diffusion. Au contraire, les structures de petite échelle en température et composition chimique n'ont pas le temps d'être diffusées. De même, le champ de vitesse offre peu de résistance au cisaillement. Ainsi la dynamo terrestre (comme celle des autres planètes) se distingue radicalement de celle du soleil et des étoiles ; la concentration de champ magnétique observée dans les taches solaires (dans un rapport de 1 à 1000) n'est possible que parce que la diffusivité thermique du plasma solaire est grande devant sa diffusivité magnétique. Au contraire, en prolongeant le champ magnétique depuis la surface de la Terre jusqu'à la surface du noyau comme il a été fait, on obtient une image assez fidèle du champ parce que, heureusement, nous savons à priori qu'il est de grande échelle. Cette comparaison entre les différentes diffusivités est centrale pour la compréhension de la dynamo terrestre, et nous donnons un dernier exemple: pour que des couches limites thermiques se développent dans l'expérience classique de chauffage d'un fluide par le bas, il est bien connu qu'il faut que la viscosité cinématique soit grande devant la diffusivité thermique. C'est le cas du manteau terrestre, par exemple. De même, dans le noyau, la forte amplitude relative de la diffusivité magnétique implique que nous devrons décrire, à la fois des couches limites visqueuses et thermiques. D'ailleurs de façon plus générale, les petites échelles du champ de vitesse (et de température) ne pourront pas être négligées (en particulier, l'apport d'énergie, qu'elle soit chimique ou thermique, a lieu à petite échelle). C'est là toute la difficulté des modélisations numériques. De plus, il faut prendre en compte des mécanismes spécifiques aux fluides en rotation rapide. En particulier, la physique des couches limites visqueuses est largement gouvernée par l'accélération de Coriolis et il se trouve que ces couches (dites d'Ekman) gardent alors un rôle important même pour des petites valeurs de la viscosité. Nous verrons que beaucoup d'études numériques se sont attachées à quantifier cet effet pour le noyau terrestre. Enfin, parce que la dimension caractéristique des structures du champ magnétique est de l'ordre du rayon du noyau, il y a aussi une composante à grande échelle des mouvements. Elle obéit à deux grands équilibres dynamiques. Il s'agit d'une part de l'équilibre, dit "magnétostrophique", entre forces magnétiques et de rotation, d'autre part de l'équilibre entre la diffusion et l'induction de champ magnétique. Il est remarquable que ces grands équilibres suffisent à expliquer les grandes propriétés du champ magnétique terrestre: la géométrie dipôle du champ correspondrait à une minimisation de la dissipation ohmique tandis que son intensité serait définie par l'équilibre entre forces magnétiques et de rotation ; quant à son orientation, on l'interprète aussi comme la conséquence du jeu entre ces deux forces [Fearn, 1994].

### 3. MODELISATION NUMERIQUE

Dans le paragraphe précédent, nous avons donné l'arrière plan des efforts de modélisation numérique de la dynamo terrestre. Une bonne partie des progrès récents s'explique uniquement par la puissance accrue des ordinateurs. Cependant, sauf pour des problèmes très idéalisés, nous ne sommes pas encore capables d'inclure dans nos modèles des échelles spatiales très distinctes. Aussi, les questions physiques les plus intéressantes (mode de transfert de l'énergie des petites aux grandes échelles et influence du champ magnétique sur les structures de petite dimension en vitesse et température) n'ont presque pas été étudiées. Le problème se simplifie cependant près du seuil. Ainsi, Cardin et Olson [1995] ont étudié le démarrage de la convection thermique, dans une sphère en rotation rapide baignée par un champ magnétique. Ils ont repris une méthode développée par
<figure><img src="image_123.png" /><figcaption>Figure 4: Représentation dans un plan équatorial des lignes de courant du fluide (a), de la distribution de température (b), et des lignes de champ magnétique induit (c). D'après Cardin et Olson [1995].</figcaption></figure>

Busse [1975] pour réduire le problème à 2 dimensions en tirant parti de la condition de Taylor-Proudman qui montre que, quand les forces de rotation sont prédominantes, les mouvements sont indépendants de la coordonnée comptée le long de l'axe de rotation. Cette approche a été récemment validée dans le cas de l'étude de la convection dans une couronne sphérique en rotation rapide, mais sans champ magnétique [Cardin et Olson, 1994] par une comparaison avec un modèle à 3 dimensions [Glatzmaier et Olson, 1993]. Dans leur dernier article, Cardin et Olson vérifient que la présence d'un champ magnétique favorise la convection, en cela qu'elle lui permet de démarrer pour un chauffage plus faible (étude linéaire) tant que l'intensité du champ imposée n'atteint pas la valeur correspondant à l'équilibre magnétostrophique entre forces magnétiques et de rotation. Pour de très faibles valeurs du champ la situation est légèrement plus compliquée et dépend du nombre de Prandtl (rapport entre viscosité et diffusivité thermique). Enfin, au démarrage de la convection, les mouvements sont organisés en rouleaux d'autant plus larges que l'intensité du champ imposée est forte. Cette conclusion se trouve renforcée par une étude de la convection développée. Comme prévu, il y a aussi des structures de petite

échelle dans la figure des mouvements et la distribution de température tandis que le champ magnétique est uniquement de grande échelle (Figure 4). Enfin, la ressemblance entre ces modèles de champ de vitesse et ceux dérivés de la variation séculaire (voir plus loin) n'a pas manquée d'être notée.

La réduction à 2 dimensions que nous venons de décrire n'est possible que quand les forces de rotation sont prédominantes. Aussi cette méthode ne peut pas être utilisée avec des champs forts. Elle ne peut pas non plus nous aider à franchir l'étape suivante qui consiste à s'affranchir du champ magnétique imposé et à modéliser un champ auto-entretenu. En effet, on s'attend à ce qu'alors la solution corresponde à un équilibre entre forces magnétiques et de rotation. On a donc cherché une autre façon de ramener le problème à deux dimensions. Dans une série d'articles publiés depuis (1964), Braginsky a proposé un modèle presque à symétrie de révolution. Un théorème ancien (1934) dû à Cowling prouve l'impossibilité de modèles de dynamos axisymétriques, mais il suffit d'une toute petite composante asymétrique pour le contrecarrer. Braginsky a étudié les ondes magnétiques planétaires se propageant dans un tel contexte (champ presque axisymétrique, vitesse presque zonale) et a montré comment paramétriser leur effet moyen (dit "effet α ") dans l'équation du terme principal à symétrie de révolution. Ainsi, le problème est à nouveau réduit à 2 dimensions.

La quasi-totalité des modèles de géodynamo dérive, à l'heure actuelle, de cette approche. Plusieurs simplifications supplémentaires sont introduites. Tout d'abord, en fixant la distribution d'anomalies de densité on écarte les difficultés associées à la faible valeur relative de la diffusivité thermique (mais on s'interdit aussi de les étudier). Ensuite le champ de vitesse $\vec{u}$ est obtenu à partir de l'équation de quantité de mouvement:

$$ \mathrm{Ro} \left( \frac{\partial}{\partial t} \vec{u} + (\vec{u} \cdot \vec{\nabla}) \vec{u} \right) + 2 (\vec{k} \wedge \vec{u}) = - \vec{\nabla} p + \vec{j} \wedge \vec{B} + \delta \rho \vec{g} + E \Delta \vec{u} \quad (2) $$

Deux petits nombres interviennent ici, le nombre d'Ekman $E = v/\Omega a^2 \approx 10^{-15}$ et le nombre de Rossby $Ro = \eta/a^2 \approx 10^{-8}$ (a rayon du noyau, $\Omega$ taux de rotation de la Terre, $\eta$ diffusivité thermique, $v$ viscosité). L'équation (2), où l'on néglige les termes en $E$ et $Ro$, correspond précisément à l'équilibre entre forces magnétiques et de rotation que nous cherchons à mettre en évidence. En négligeant ces termes, on écrirait une équation dite "diagnostique" pour le champ de vitesse que l'on pourrait calculer pour tout champ magnétique B et distribution d'anomalies de densité $\delta p$. En reportant cette expression pour le champ de vitesse dans l'équation d'induction, nous n'aurions plus qu'à suivre l'évolution dans le temps d'une seule inconnue: le champ magnétique. Nous serions ramenés à l'étude de la seule équation de l'induction, c'est à dire de celle qui présente le moins de difficultés en raison de la forte diffusivité magnétique.

C'est pourquoi, le travail de J.B. Taylor [1963], qui étudia sous quelle condition il est possible de négliger les termes en $E$ et $Ro$, prît tant d'importance. Il montra que l'équation (2) réduite n'a de solution qu'elle si la composante azimutale des forces magnétiques $(\vec{j} \wedge \vec{B}) \cdot e_\phi$. $(e_\phi$ vecteur unitaire le long d'un parallèle) est de moyenne nulle sur toute surface cylindrique centrée sur l'axe de rotation de la Terre. Sinon, elle accélère la rotation de ces cylindres. Il y alors contradiction avec l'hypothèse de départ. Il faut soit réintroduire le terme en $Ro$, soit tenir compte de la modification de la condition aux limites pour le champ de vitesse dans l'intérieur du noyau due à la couche-limite visqueuse attachée au manteau, qui joue un rôle beaucoup plus important
<figure><img src="image_124.png" /><figcaption>Figure 5: Description de la nature des solutions pour un modèle αω particulier (intensité du champ poloidal $B_p$ en fonction du nombre dynamo D représentant l'amplitude du forçage). Au seuil (D et $B_p$ faibles, branche a), l'amplitude du champ est limitée par la viscosité. Sur une seconde branche (b), il y a équilibre entre forces magnétiques et de rotation.</figcaption></figure>

que les forces visqueuses de volume parce qu'elle intervient à l'ordre $E^{1/2}$ et non pas E. La plupart des auteurs n'ont traité que cette dernière possibilité: les forces magnétiques de volume sont compensées par un frottement visqueux au bord. Dans ce cadre, plusieurs types de solutions deviennent alors possibles. Au démarrage de l'action dynamo, l'amplitude du champ magnétique est toujours d'ordre $E^{1/4}$ (régime dit d'Ekman). Pour des forçages plus grands, on peut rechercher des solutions telles que la géométrie du champ satisfasse la condition écrite par Taylor (d'où le nom de régime de Taylor). Enfin, on peut imaginer des solutions mixtes comme celle proposée par Braginsky [1975]: en s'alignant le long de l'axe de rotation de la Terre, le champ magnétique diminue l'amplitude des forces magnétiques zonales de moyenne non nulle sur les cylindres géostrophiques mais l'intensité du champ, dans son modèle, dépend quand même de la viscosité. Ces différentes solutions se distinguent par leur comportement asymptotique dans la limite des très petites visco sités, limite pertinente pour le noyau de la Terre. Il n'était pas possible de le mettre en évidence numériquement tant que l'on n'était pas capable de travailler à de très faibles viscosités. Or Jault [1995] vient de montrer que réintroduire un terme Ro (avec $Ro \ll E^{1/2}$) le permet.

Un régime de Taylor (c'est à dire l'équilibre magnétostrophique) a pu ainsi être mis en évidence pour tous les modèles étudiés, y compris celui qui avait été proposé pour illustrer la théorie de Braginsky. Il semble donc que ce soit la solution générale pour les faibles viscosités et un apport d'énergie suffisamment grand (Figure 5). Enfin, la même solution peut être retrouvée en négligeant complètement tout terme visqueux mais en conservant le terme d'accélération en Ro, qui permet d'atteindre l'état final pourtant indépendant de Ro. Soulignons que la réussite de cette méthode vient de ce qu'elle est calquée sur le mécanisme qui permet au champ magnétique terrestre d'atteindre son équilibre. Tout écart à un état défini par la condition de Taylor (ou une condition modifiée pour tenir compte de la non-sphéricité de la frontière noyau-manteau) provoque des ondes de torsion d'Alfven qui permettent au champ d'adopter une nouvelle géométrie d'équilibre. Jault et Le Mouël [1993] ont d'ailleurs proposé que ces mouvements (rotation de cylindres centrés sur l'axe de rotation de la Terre) soient respon-

sables des variations du moment cinétique du noyau et donc par contrecoup de la Terre solide. Nous verrons que cette idée s'est avérée fructueuse.

## 4. LES MOUVEMENTS A LA SURFACE DU NOYAU RESPONSABLES DE LA VARIATION SECULAIRE

### 4.1. Hypothèse du flux gelé

Disposant, comme nous l'avons vu, d'une description $\{g_n^m(t), h_n^m(t)\}$ du champ magnétique principal et de sa variation séculaire variant sur près de trois siècles à la base du manteau, et partant de ce que cette variation séculaire semble bien être la signature des mouvements associés au mécanisme dynamo (voir ci-dessus), nous allons maintenant voir comment il est possible d'obtenir des informations sur la trace de ces mouvements à la surface du noyau.

Les variations du champ magnétique dans le noyau liquide et conducteur sont gouvernées par l'équation de l'induction. Explicitant la composante radiale de cette équation au sommet du noyau (considéré comme sphérique) et utilisant la condition aux limites $u_r = 0$, nous obtenons:

$$ \frac{\partial B_r}{\partial t} = - \vec{\nabla}_H \cdot (\vec{u} B_r) + \frac{\eta_m}{r} \Delta (r B_r) \quad (3) $$

où $\vec{\nabla}_H$ est la composante horizontale de $\vec{\nabla}$ et $\eta_m = (\mu\sigma)^{-1}$ est la diffusivité magnétique ($\sigma$ étant la conductivité du noyau). Oubliant un instant le terme de diffusion (dernier terme), on constate que le champ magnétique n'intervient dans (3) que par le biais de $B_r$, de ses dérivées horizontales et de sa dérivée temporelle $(\delta B_r / \delta t)$, toutes quantités connues (dans la limite de précision des modèles géomagnétiques) par continuité lors du passage de la base du manteau à la surface du noyau. Le terme de diffusion, qui fait par contre également intervenir des dérivées radiales de $B_r$ inaccessibles à l'observation, est généralement considéré comme négligeable. On peut avoir recours à deux types de raisonnement pour justifier de cette hypothèse (dite du "flux gelé").

On peut d'abord remarquer que si l'on néglige la diffusion dans (3), on doit alors avoir conservation au cours du temps de la quantité intégrale (dite "flux magnétique") $\Phi = \iint B_r dS$ sur toute surface délimitée par une courbe $\Gamma_{B_r=0}$ au sommet du noyau. Nous avons vu qu'il était possible de suivre l'évolution de ce type de surfaces et donc d'apprécier l'évolution des flux $\Phi$, au moins d'un point de vue global (l'identification de chacune des taches n'étant pas toujours aisée). Bloxham et al. [1989] ont constaté que ces flux $\Phi$ ne variaient pas de manière significative sauf, semble-t-il, pour une "tache" à l'aplomb du sud de l'Afrique (Figure 1). L'estimation de la variation constatée est cependant contestée [Backus, 1988] car elle est à la limite de la précision que l'on peut accorder aux modèles géomagnétiques. Dans le doute, on peut alors tenter d'apprécier l'ordre de grandeur des dérivées radiales de $B_r$ à partir de la physique du problème et estimer directement l'influence du terme diffusif (les dérivées horizontales observées de $B_r$ ne donnant lieu à aucune contribution significative). Bloxham [1986] a par exemple proposé un mécanisme d'expulsion de flux pouvant apparemment rendre compte de la variation de flux décrite précédemment. Cependant son modèle "plan" n'est pas compatible avec la situation réelle pour des raisons géométriques [Jault, 1990]. En fait il semble difficile de rendre compte d'une telle observation. Ainsi, si l'on se
place dans l'hypothèse réaliste (voir Cardin et Olson [1995], paragraphe précédent) d'un mouvement de grande échelle (L ~ 10⁶ m), le terme diffusif en surface est essentiellement dû à un "effet de peau" imposé par la présence du manteau isolant [e.g., Jault et Le Mouël, 1991a]. L'ordre de grandeur du terme diffusif relativement à la variation séculaire (δBᵣ) / (δt) étant alors de 0 (√ηTₛV/L) ~ 0.1, si l'on utilise la valeur σ = 7.10⁵ (Ωm)⁻¹ particulièrement bien contrainte pour la conductivité [e.g., Poirier, 1994], et la valeur TₛV ~ 200 ans que l'on peut déduire de l'étude statistique de Hulot et Le Mouël [1994], précédemment décrite. La diffusion est donc d'ampleur modeste et peut être légitimement négligée. Par contre il importe de signaler que, du fait même «effet de peau», le terme diffusif n'est pas tout à fait négligeable lorsque l'on considère la composante horizontale de l'équation de l'induction (complémentaire de (3)) [Jault et Le Mouël, 1991a]. C'est pourquoi nous ne l'emploierons pas. Enfin, on peut imaginer que le mouvement soit très fortement stratifié, mais même dans cette situation, il est difficile de ne pas retomber sur une propriété de type «flux gelé» [Braginsky et Le Mouël, 1993].

## 4.2. Calcul de mouvements tangentiellement géostrophiques

Nous plaçant désormais dans l'hypothèse d'un mouvement de grande échelle, négligeant la diffusion, connaissant tous les termes autres que le mouvement ū dans l'équation (3), que nous baptiserons alors "équation de la variation séculaire", nous sommes presque à même de déterminer les mouvements à la surface du noyau. Mais la solution au problème mathématique ainsi posé n'est pas unique. Heureusement, une information supplémentaire sur la nature des mouvements peut être obtenue à partir de l'équation (2) qui, comme il a déjà été dit, peut se réduire à l'équation diagnostique (écrite maintenant sous forme dimensionnée):

$$2\rho_0 (\vec{\Omega} \wedge \hat{u}) = -\vec{V}p + \vec{j} \wedge \vec{B} + \delta\rho \vec{g} \quad (4)$$

Nous faisons donc l'hypothèse, cohérente avec les résultats obtenus par Jault [1995, voir ci-dessus], que les mouvements correspondent à un régime proche de celui de Taylor. Considérée à la surface du noyau, l'équation (4) peut encore être un peu simplifiée. Utilisant toujours la même approche en couche limite diffusive de Jault et Le Mouël [1991a], et pour les valeurs du champ magnétique observé, on peut en effet montrer que la force de Lorentz j ∧ B est d'environ 10⁻³ la force de Coriolis. Sachant enfin que la composante non radiale de δρ⃗g est elle aussi négligeable [Hulot et al., 1995], on peut alors écrire la composante horizontale de (4) sous la forme:

$$2\rho_0 (\vec{\Omega} \wedge \hat{u})_H = -\vec{V}_{HP} \quad (5)$$

ce qui constitue l'hypothèse tangentiellement géostrophique [Le Mouël, 1984]. Armée de cette contrainte, la solution de l'équation de la variation séculaire est alors presque unique, l'indétermination ne portant plus que sur des mouvements circulant dans des régions limitées de la surface du noyau [Backus et Le Mouël, 1986]. Il devient possible de déterminer de manière unique un mouvement de grande échelle à partir des modèles de champ et de variation séculaire pour une époque donnée [Gire et Le Mouël, 1990].

Un tel mouvement est représenté par deux champs scalaires S et T, eux-mêmes développés en harmoniques sphériques, et tels

<figure><img src="image_125.png" /><figcaption>Figure 6: Vue polaire (Nord) du scalaire toroidal T (partie non zonale), associé au mouvement symétrique rendant compte des observations magnétique pour 1980. Le mouvement est symétrique par rapport à l'équateur (non visible sur cette représentation) et possède un plan de symétrie à environ 90° de Greenwhich (donc quasi vertical sur cette représentation). On notera la parenté avec la Figure 4a. D'après Hulot et al., [1990].</figcaption></figure>

que ū = c⃗V_HS - c⃗n ∧ ⃗V_HT, où c est le rayon du noyau et n le vecteur unitaire radial. Le scalaire poloidal S décrit la partie du mouvement qui implique des échanges de fluide entre la surface du noyau et les profondeurs (les écoulements se faisant perpendiculairement aux lignes d'isovaleurs de S, le fluide converge et plonge là où S est maximum, et émerge et diverge là où au contraire S est minimum). Le scalaire toroidal T décrit des circulations de fluide restant en surface (les écoulements se faisant parallèlement aux lignes d'isovaleurs de T, dans le sens trigonométrique autour d'un maximum, inverse autour d'un minimum). Si cette écriture permet de relativement bien se représenter le mouvement, l'hypothèse tangentiellement géostrophique impose cependant des liens entre S et T. En effet à partir de l'équation (5) on montre que l'on a aussi ū = (2ρ₀Ωc cosθ)⁻¹n ∧ ⃗V_Hp. Ainsi le mouvement ne dépend en fait que d'un seul champ scalaire, le champ de surpression. Ceci amène du reste un résultat qui nous sera fort utile pour la suite: le champ de surpression à la base du manteau est calculé en même temps que le mouvement proprement dit. Signalons pour finir que, compte tenu de notre connaissance limitée au degré 13 pour le champ principal et au degré 7 pour la variation séculaire, le mouvement et la pression ne peuvent en aucun cas être calculés au delà du degré 4 [Hulot et al., 1992].

## 4.3. Non-stationnarité des mouvements et secousses géomagnétiques

Les premiers calculs de mouvements rigoureusement géostrophiques ont été réalisés pour des époques voisines de 1980 [Gire et Le Mouël, 1990; Hulot et al., 1990] et ont immédiatement fait apparaître une tendance des mouvements à présenter de fortes symétries: une symétrie équatoriale doublée d'une symétrie centrale se traduisant par la présence d'un plan méridional de symétrie (Figure 6). Les calculs montraient qu'il était même possible d'expliquer les données magnétiques par un mouvement présentant exactement ces symétries et il était tentant de voir là
<figure><img src="image_126.png" /><figcaption>Figure 7: Variations décennales de la durée du jour observées (carrés vides), prédits par les mouvements de Hulot et al. [1993] (carrés pleins). Unités en ms.</figcaption></figure>

une éventuelle propriété intrinsèque de la géodynamo. Si des calculs ultérieurs ont dans un premier temps semblé confirmer cette impression (pour l'époque 1900, *Hulot* [1922]), on sait maintenant grâce à une étude portant sur 150 ans d'observations que les choses ne sont pas aussi simples et, en particulier, que les mouvements sont susceptibles d'évoluer énergétiquement et morphologiquement beaucoup plus rapidement que ce que l'on imaginait [*Hulot et al.*, 1994]. Il faut donc concevoir les mouvements comme étant fondamentalement non stationnaires à l'échelle décennale.

Cette non-stationnarité semble d'ailleurs liée à un aspect de la variation séculaire que nous avons déjà décrit et que avons pour l'instant ignoré: l'occurrence de secousses géomagnétiques. Celles-ci semblent en effet associées aux changements observés dans le régime des mouvements à la surface du noyau [*Hulot et al.*, 1994]. Ce résultat est d'autant plus intéressant qu'une seconde étude de nature plus théorique [*Jault et Le Mouël*, 1994] a permis de montrer que seuls des mouvements intéressant au moins une couche d'environ 100 km d'épaisseur étaient susceptibles de générer de telles secousses. Ainsi, si la nature exacte de ces secousses reste encore mal comprise, il apparaît comme probable qu'elles trouvent leur origine dans les mêmes mouvements que ceux qui génèrent la variation séculaire proprement dite.

## 4.4. *Non-stationnarité des mouvements et rotation de la Terre*

Le fait que les mouvements à la surface du noyau ne soient pas stationnaires suggère que l'ensemble des mouvements animant le

noyau, et en particulier ceux contribuant au moment cinétique nucléaire, soient également non stationnaires. Une conséquence mécanique immédiate d'une telle éventualité serait l'apparition, à l'échelle décennale, de perturbations dans la rotation de la Terre (terme du reste mal nommé, car il s'agit en fait de la rotation du manteau, les mesures astronomiques étant faites relativement à celui-ci). On peut tenter de prédire ces perturbations à partir des mouvements calculés précédemment. Deux grands types de calculs peuvent alors être menés, soit que l'on s'appuie sur le principe de conservation du moment cinétique, ce qui dispense de connaître le mécanisme exact par lequel le noyau est couplé au manteau mais requiert une hypothèse sur la façon dont les mouvements en surface peuvent être extrapolés en profondeur, soit que l'on propose un mécanisme de couplage, ce qui permet de n'avoir à faire dans la plupart des cas que peu d'hypothèse sur ce qu'il se passe dans le noyau profond. Par ailleurs, le problème de la rotation se découle naturellement en un problème dit "axial", qui intéresse les irrégularités de la vitesse de rotation (et donc de la durée du jour), et un problème dit "équatorial", qui concerne les variations de l'orientation de l'axe de rotation du manteau (et donc de la polhodie, mouvement de la trace du pôle de rotation sur le sol). On peut évidemment aborder chacun de ces deux problèmes par l'une et l'autre approches.

Le premier résultat à avoir ainsi été obtenu a déjà été brièvement évoqué à la fin du paragraphe concernant la modélisation numérique. Partant de l'hypothèse que la convection dans le noyau est proche du régime de Taylor, nous avons vu en effet que tout léger départ à ce régime se traduisait par l'apparition de mouvements accélérés en cylindres centrés sur l'axe de rotation de la Terre. Ceci se produit en particulier s'il existe un couplage noyau-manteau et ce, quel que soit le mécanisme de couplage [*Jault et al.*, 1988; *Jault et Le Mouël*, 1989, 1993]. En prolongeant dans l'ensemble du noyau la trace, visible en surface, de ces mouvements en cylindre, on peut alors estimer le moment cinétique et prédire les variations décennales de la durée du jour. Ce calcul est d'autant plus facile à mener que seuls interviennent les composantes toroidales zonales de degré 1 et 3 des mouvements calculés à partir des données magnétiques (les autres composantes n'étant pas porteuses de moment cinétique). La comparaison entre prédiction et observation est satisfaisante même si l'accord n'est pas parfait pour les époques les plus anciennes (Figure 7).

On peut alors tenter d'aborder le problème axial en termes de couplage noyau-manteau. Ce dernier peut a priori être visqueux, électromagnétique, topographique, ou gravitationnel. Il est hors de propos de reprendre la discussion portant sur chacun de ces couplages, et nous nous bornerons à résumer l'opinion qui prévaut aujourd'hui [e.g., *Jault et Le Mouël*, 1991b, 1993]. L'extrêmement faible viscosité [e.g., *Poirier*, 1994] permet d'éliminer le couplage visqueux. Le couplage électromagnétique semble plutôt trop faible. Par contre le couplage topographique qui correspond à l'action de la surpression sur la topographie de la frontière noyau-manteau, semble beaucoup trop grand (de un à deux ordres de grandeur, [*Hinderer et al.*, 1990; *Jault et Le Mouël*, 1990)) si l'on se fie aux topographies couramment admises [e.g., *Morelli et Dziewonski*, 1987] et aux pressions calculées avec les mouvements. Enfin, le couplage gravitationnel demeure mal connu pour la raison que les hétérogénéités à l'origine de la convection dans le noyau sont difficilement estimables. C'est pourquoi *Jault et Le Mouël* [1989] ont imaginé soit que le couple gravitationnel compense la majorité du couple topographique, soit que ce couple topographique puisse être plus faible que ce que donne un calcul nécessairement approché. En effet, en modifiant légèrement (dans les barres d'erreur) la position relative des mouvements et de la topographie, on peut diminuer la valeur de ce couple du facteur requis. Cette astuce
mathématique correspondrait en outre à la condition physique que les mouvements se maintiennent "orthogonaux" à la topographie [*Jault et Le Mouël*, 1990]. C'est ainsi que le noyau se maintiendrait proche de la condition de Taylor en présence d'une topographie.

Le problème équatorial a fait jusqu'à présent l'objet de relativement moins d'attention. Une raison à cela vient de ce qu'il n'existe pas, à l'échelle décennale, de réaction du noyau d'ampleur comparable à celle, en anneaux, que l'on observe lorsqu'apparaît un couple axial. Alors que dans le cas axial le noyau se voit contraint de modifier son moment cinétique par une accélération relative au manteau, il peut se contenter, dans le cas équatorial, de réagir en modifiant le couple non nul des forces de Coriolis. Or ceci ne nécessite qu'une faible modification quasi-stationnaire des mouvements. En d'autres termes, nous ne disposons pas d'une contrainte pour le problème équatorial comparable à la contrainte de Taylor pour le problème axial. Il est donc difficile de mener un raisonnement analogue, de type conservation du moment cinétique. Par contre les calculs de couples d'interaction noyau-manteau peuvent être menés dans le cas présent, sans difficulté supplémentaire par rapport au cas axial. De nouveau le couplage visqueux est beaucoup trop faible et le couplage électromagnétique pas assez fort [*Greff-Leffitz et Legros*, 1995], de sorte que le candidat privilégié reste à priori le couplage topographique [*Hinderer et al.*, 1987]. Curieusement les premiers calculs menés avec les mêmes topographies que précédemment donnent un couple un peu trop faible [*Greff-Leffitz et Legros*, 1995; *Hulot et al.*, 1995). Toutefois l'ampleur exacte des variations décennales de la polhodie, autrement dit du "Markovitz wobble", reste relativement mal connue et il n'est pas exclu qu'en tenant compte des composantes de plus petite échelle et donc inaccessibles à l'observation, de la topographie et des mouvements à la surface du noyau, on puisse faire s'accorder le couple topographique et l'observation. Ceci suppose néanmoins que l'on ne fasse jouer aucun rôle au couplage gravitationnel. Or celui-ci, comme dans le cas axial, reste l'un des grands paramètres inconnus du problème et il est incontestable qu'il joue un rôle important au moins dans le cas simplifié d'une topographie elliptique [*Hulot et al.*, 1995].

Ainsi, que l'on considère le problème axial ou le problème équatorial, il apparaît que les couplages les plus importants sont les couplages topographique et gravitationnel et, dans une moindre mesure, le couplage électromagnétique. Il est clairement nécessaire en outre de tenir compte de ces couplages simultanément. Il reste enfin que les conclusions que l'on peut tirer de l'étude du problème axial d'une part, du problème équatorial d'autre part, se doivent évidemment d'être cohérentes. C'est en considérant les deux problèmes simultanément que de nouveau progrès pourront donc être réalisés.

## REFERENCES

Achache, J., V. CoOurtillot, J. Ducruix and J. L. Le Mouël, The late 1960's secular variation impulse: further constraints on deep mantle conductivity, *Phys. Earth Planet. Inter.*, 23:72-75, 1980.

Alexandrescu, M., D. Gibert, G. Hulot, J. L. Le Mouël, G. Saracco, Detection of geomagnetic jerks using wavelet analysis, *J. Geophys. Res.*, sous presse, 1995.

Alldredge, L. R., A discussion of impulses and jerks in the geomagnetic field, *J. Geophys. Res.*, 89:4403-4412, 1984.

Backus, G. E., Application of mantle filter theory to the magnetic jerk of 1969, *Geophys. J. R. Astr. Soc.*, 74:713-746, 1983.

Backus, G. E., Bayesian inference in geomagnetism, *Geophys. J.*, 92:125-142, 1988.

Backus, G. E., J. L. Le Mouël, The region on the core-mantle boundary where a geostrophic velocity field can be determined from frozen-flux magnetic data, *Geophys. J. R. Astr. Soc.*, 85:617-628, 1986.

Backus, G. E., R. H. Estes, D. Chinn, R. A. Langel, Comparing the jerk with other global models of the geomagnetic field from 1960 to

1978, *J. Geophys. Res.*, 92:3615-3622, 1987.

Bloхjam, J., The expulsion of magnetic flux from the Earth's core, *Geophys. J. R. Astr. Soc.*, 87:669-678, 1986.

Bloхjam, J., A. Jackson, Time-dependent mapping of the magnetic field at the core-mantle boundary, *J. Geophys. Res.*, 97:19537-19564, 1992.

Bloхjam, J., D. Gubbins, A. Jackson, Geomagnetic secular variation, *Phil. Trans. R. Soc. Lond.*, A329:415-502, 1989.

Braginsky, S. I., Self-excitation of a magnetic fluid during the motion of a highly conducting fluid, *Soviet Phys. JETP*, 20:726-735, 1964.

Braginsky, S. I., Nearly axially symmetric model of the hydromagnetic dynamo of the Earth I, *Geomagnetism and Aeronomy*, 15:122-128, 1975.

Braginsky, S. I., J. L. Le Mouël, Two-scale model of a geomagnetic field variation, *Geophys. J. Int.*, 112:147-158, 1993.

Bucur, I., The direction of the terrestrial magnetic field in France during the last 21 centuries, recent progress, *Phys. Earth Planet. Inter.*, 87:95-109, 1994.

Busse, F. H., A model of the Geodynamo, *Geophys. J. Roy. Astr. Soc.*, 42:437-459 1975.

Cardin, P., P. Olson, Chaotic thermal convection in a rapidly rotating spherical shell: Consequences for flow in the Outer core, *Phys. Earth Planet. Interiors*, 82:235-259, 1994.

Cardin, P., P. Olson, The influence of toroidal magnetic field on thermal convection in the core, *soumis à Earth Planet. Sciences Letters*, 1995.

Cohen, Y., Traitements et interprétations de données spatiales en géomagnétisme: étude des variations latérales d'aimentation de la lithosphère terrestre, *Thèse, Université Paris VII*, 1989.

Cohen, Y., J. Achache, New global vector magnetic anomaly maps derived from magsat data, *J. Geophys. Res.*, 95:10783-10800, 1990.

Counil, J. L., Y. Cohen, J. Achache, The global continent-ocean magnetization contrast, *Earth Planet. Sci. Lett.*, 103:354-364, 1991.

Courtillot, V., J. L. Le Mouël, Time variations of the Earth's magnetic field: from daily to secular, *Ann. Rev. Earth Planet. Sci.*, 16:389-476, 1988.

Courtillot, V., J. Ducruix, J. L. Le Mouël, Sur une accélération récente de la variation séculaire du champ magnétique terrestre, *C. R. Acad. Sci. Paris*, 287:1095-1098, 1978.

Courtillot, V., J. P. Valet, G. Hulot, J. L. Le Mouël, The Earth's magnetic field: which geometry?, *Eos Trans. AGU*, 73:337;340-342, 1992.

Cowling, T. G., The magnetic field of sunspots, *Mon. Not. R. Astr. Soc.*, 94:39-48, 1934.

Daly, L., M. Le Goff, An updated and homogeneous world secular variation database. 1-Smoothing of the Archeomagnetic results, *soumis à Phys. Earth Planet. Inter.*, 1995.

Fearn, D. R., Non linear planetary dynamos, in Solar and Planetary dynamos, Vol. 2, édité par Proctor, Matthews et Rucklidge, Cambridge University Press, 1994.

Gire, C., Sur la variation séculaire du champ magnétique terrestre et les mouvements des couches externes du noyau fluide, *Thèse, Université Paris VII*, 1985.

Gire, C., J. L. Le Mouël, Tangentially geostrophic flow at the core-mantle boundary compatible with the observed geomagnetic secular variation: the large-scale component of the flow, *Phys. Earth Planet. Inter.*, 59:259-287, 1990.

Glatzmaier G. Olson P. "Highly supercritical thermal convection in a rotating spherical shell: centrifugal vs radial gravity" *Geophysical and astrophysical fluid dynamics* 70:113-126, 1993.

Greff-Leffitz, M., H. Legros, Core-mantle coupling and polar motion, *soumis à Phys. Earth Planet. Inter.*, 1995.

Gubbins, D., P. Kelly, Persistent patterns in the geomagnetic field over the past 2.5 Myr, *Nature*, 365:829-832, 1994.

Hinderer, J., C. Gire, H. Legros, J. L. Le Mouël, Geomagnetic secular variation, core motions and implications for the Earth's wobble, *Phys. Earth Planet. Inter.*, 49:121-132, 1987.

Hinderer, J., D. Jault, H. Legros, J. L. Le Mouël, Core-mantle topographic torque: a spherical harmonic approach and implications for the excitation of the Earth's rotation by core motions, *Phys. Earth Planet. Inter.*, 59:329-341, 1990.

Hulot, G., Observations géomagnétiques et Géodynamo, *Thèse, Université Paris VII*, 1992.

Hulot, G., J. L. Le Mouël, A statistical approach to the Earth's main magnetic field, *Phys. Earth Planet. Inter.*, 82:167-183, 1994.

Hulot, G., J. L. Le Mouël, D. Jault, The flow at the core-mantle boundary: symmetry properties, *J. Geomag. Geoelectr.*, 42:857-874, 1990.

Hulot, G., J. L. Le Mouël, J. Wahr, Taking into account truncation problems and geomagnetic model accuracy in assessing computed flows at the core-mantle boundary, *Geophys. J. Int.*, 108:224-246, 1992.
Hulot, G., M. Le Huy, J. L. Le Mouël, Secousses (jerks) de la variation séculaire et mouvements dans le noyau terrestre, *C. R. Acad. Sci. Paris, Série II*, 317:333-341, 1993.

Hulot, G., L. Hongre, A. Khokhlov, J. L. Le Mouël, From secular to short Paleomagnetic time scales: the Earth's main magnetic field as a stochastic process, *Abstract book, Fourth SEDI symposium*, 1994.

Hulot, G., M. Le Huy, J. L. Le Mouël, Influence of core flows on the decade variations of the polar motion, *soumis à Geophys. Astrophys. Fluid Dyn.*, 1995.

Jault, D., Variation séculaire du champ géomagnétique et fluctuations de la longueur du jour, *Thèse, Université Paris VII*, 1990.

Jault, D., Model Z by computation and the Taylor's condition, *Geophys. Astrophys. Fluid Dyn.*, sous presse, 1995.

Jault, D., J. L. Le Mouël, The topographic torque associated with a tangentially geostrophic motion at the core surface and inferences on the flow inside the core, *Geophys. Astrophys. Fluid Dyn.*, 48:273-296, 1989.

Jault, D., J. L. Le Mouël, Core-mantle boundary shape: constraints inferred from the pressure torque acting between the core and the mantle, *Geophys. J. Int.*, 101:233-241, 1990.

Jault, D., J. L. Le Mouël, Physical properties at the top of the Earth's core and core surface motions, *Phys. Earth Planet. Inter.*, 68:76-84, 1991a.

Jault, D., J. L. Le Mouël, Exchange of angular momentum between the core and the mantle, *J. Geomag. Geoelectr.*, 43:111-129, 1991b.

Jault, D., J. L. Le Mouël, Circulation in the liquid core and coupling with the mantle, *Adv. Space Res.*, 13:(11)221-(11)233, 1993.

Jault, D., J. L. Le Mouël, Does secular variation involve motions in the deep core ?, *Phys. Earth Planet. Inter.*, 82:185-193, 1994.

Jault, D., C. Gire, J. L. Le Mouël, Westward drift, core motions and exchanges of angular momentum between core and mantle, *Nature*, 333:353-356, 1988.

Le Mouël, J. L., Outer-core geostrophic flow and secular variation of the Earth's geomagnetic field, *Nature*, 311:734, 1984.

Le Mouël, J. L., J. P. Poirier, H. Assoumani, Sur le refroidissement du noyau terrestre, *C. R. Acad. Sci. Paris, Série II*, 319:1497-1501, 1994.

Prévot, M., M. Perrin, Intensity of the Earth's magnetic field since Precambrian from Thellier-type paleointensity data and inferences on the thermal history of the core, *Geophys. J. Int*, 108:613-620, 1994.

McLeod, M. G., Signals and noise in magnetic observatory annual means: mantle conductivity and jerks, *J. Geophys. Res.*, 97:17261-17290, 1992.

Morelli, A., A. M. Dziewonski, Topography of the core-mantle boundary and lateral homogeneity of the liquid core, *Nature*, 325:678-683, 1987.

Poirier, J. P., Physical properties of the Earth's core, *C. R. Acad. Sci. Paris, Série II*, 318:341-350, 1994.

Poirier, J. P., J. L. Le Mouël, Does infiltration of core material into the lower mantle affect the observed geomagnetic field?, *Phys.Earth Planet.Inter.*, 1992.

Quidelleur, X., J. P. Valet, V. Courtillot, G. Hulot, Long-term geometry of the geomagnetic field for the last five million years: an updated secular variation database, *Geophys. Res. Lett.*, 21:1639-1642, 1994.

Stewart, D. N., Geomagnetic impulses and the electrical conductivity of the lower mantle, *Ph. D. Thesis, The University of Leeds, Department of Earth Sciences*, 1991.

Taylor, J. B., The magneto-hydrodynamics of a rotating fluid and the Earth's dynamo problem, *Proc. R. Soc. Lond.*, A274:274-283, 1963.
# LES EMISSIONS GAZEUSES EN RÉGION TROPICALE

Robert DELMAS, Jean Pierre LACAUX et Bernard CROS

Laboratoire d'Aérologie (URA CNRS 354), Université Paul Sabatier,
118 Route de Narbonne, 31 062 Toulouse Cedex

Hélène CACHIER

Centre des Faibles Radioactivités, Domaine du CNRS, 91 198 Gif sur Yvette Cedex

L'influence des zones tropicales sur la composition chimique de l'atmosphère globale, peu étudiée jusqu'à la fin des années 70, est désormais considérée comme un des paramètres les plus importants pour l'étude des changements globaux. Cet article présente un panorama des échanges biosphère-atmosphère dans la zone intertropicale et de leurs conséquences sur la chimie atmosphérique. Il est illustré par des exemples s'appuyant sur les travaux réalisés notamment en Afrique, par des chercheurs français. Les sources de composés gazeux dans les régions tropicales sont principalement la combustion de biomasse et les émissions biogéniques des écosystèmes terrestres. Les composés à l'état de trace produits par la biosphère continentale influencent la composition chimique et le pouvoir oxydant de l'atmosphère, à l'échelle régionale et globale, au travers de processus homogènes et hétérogènes. L'impact radiatif des émissions gazeuses et particulaires est complexe mais il est clairement établi que le déstockage de carbone de la biosphère terrestre lors de la destruction des forêts tropicales contribue à l'effet de serre additionnel.

The influence of tropical regions on the chemical composition of the atmosphere is now considered as one of the most important parameter of global change studies. This paper presents an overview of biosphere-atmosphere exchanges in the tropics and of their consequences on atmospheric chemistry. It is illustrated by examples taken from research activities conducted in Africa by French scientists. Biomass burning and biogenic emissions from terrestrial ecosystems are the main sources of trace gases and aerosols in the tropics. Tace compounds produced by the continental biosphere affect the chemical composition and the oxidation capacity of the atmosphere, at regional and global scale, through homogeneous and heterogeneous processes. The radiative impact of gas and particle emissions is complex however it is clear that the carbon destorage from the terrestrial biosphere, due to deforestation, enhances the additional greenhouse effect.

## 1. INTRODUCTION

La composition de l'atmosphère globale est d'abord le résultat de l'activité de la biosphère. On sait depuis longtemps qu'il en est ainsi pour les principaux constituants atmosphériques: l'oxygène moléculaire est un "produit secondaire" de la photosynthèse des plantes et l'azote est produit par les bactéries des sols pendant la dénitrification. Sans la régénération continue de l'oxygène et de l'azote moléculaire par la biosphère, l'atmosphère irait vers un équilibre thermodynamique représenté par le dioxyde de carbone et les nitrates. Les échanges gazeux avec la biosphère ont donc une influence prépondérante sur les propriétés physiques et chimiques de l'atmosphère de la Terre.

Les écosystèmes tropicaux présentent la biogéochimie la plus dynamique mais la moins connue. Les tropiques représentent 40% de la biomasse globale et environ 60% de la production primaire nette annuelle et globale (PPN). Cette productivité est caractérisée par une diversité des espèces végétales sans équivalent dans aucune autre région. Elle est alimentée par d'importants apports d'énergie solaire et par la transformation de cette énergie en chaleur latente et sensible, qui recycle l'humidité en régimes météorologiques complexes mais prévisibles. Ces facteurs déterminent également le rôle important des régions intertropicales sur la dynamique atmosphérique globale.

La chaleur et l'humidité commandent la plupart des cycles biogéochimiques y compris la production primaire nette et le cycle du carbone. Les flux de gaz traces leur sont intimement associés, ils sont donc souvent tributaires des variations saisonnières des paramètres climatiques. Ces émissions incluent des composés qui peuvent affecter le bilan radiatif global comme le méthane [CH₄] et l'oxyde nitreux [N₂O] ou le bilan régional des oxydants comme les hydrocarbures réactifs, les oxydes d'azote, le monoxyde de carbone [CO]. Par ailleurs les produits dérivés des transformations photochimiques des constituants primaires

comme le nitrate de peroxyacétyl [PAN], l'ozone [O₃] et les espèces radicalaires associées affectent le temps de vie des gaz radiativement actifs. Enfin, ces oxydants déterminent la formation de composés acides qui par leur dépôt, sec ou humide, peuvent provoquer des effets physiologiques rétroactifs sur les émissions (O₃, pluies acides).

Les campagnes conduites au cours des quinze dernières années ont confirmé l'importance des tropiques sur la composition globale et régionale et la chimie de l'atmosphère. Les campagnes Quemadas au Brésil en 1979 et 1982 ont montré l'importance des feux de biomasse. Les missions GTE/ABLE 2A et 2B (Global Tropospheric Experiments) en Amazonie ont défini plus tard la composition chimique régionale du bassin de l'Amazone pendant la saison humide et la saison sèche. En Afrique, des études ont été réalisées au Congo (1988) et en Côte d'Ivoire (1991) dans le cadre du programme DECAFE (Dynamique Et Chimie de l'Atmosphère en Forêt Equatoriale), études portants principalement sur les émissions naturelles ou par les feux de constituants mineurs et la formation d'oxydants au niveau des écosystèmes tropicaux africains. Récemment des études supplémentaires ont été conduites en Amérique du Sud pendant les missions STARE/TRACE et en Afrique du Sud dans le cadre de TRACE/SAFARI (1992).

## 2. LES SOURCES DE CONSTITUANTS MINEURS ATMOSPHÉRIQUES DANS LES RÉGIONS TROPICALES

Bien que la pollution anthropique, liée à la combustion d'énergie fossiles et aux activités humaines en général, se soit considérablement accrue dans des régions à hautes densités de population (Asie, Amérique du Sud), dans la plupart des zones tropicales, les principales sources de gaz trace et d'aérosols résultent des échanges biosphère-atmosphère. Elles incluent les émissions naturelles des sols et de la végétation dans des
<table><thead><tr><th rowspan="2">Espèces</th><th rowspan="2">taux d'émission (10<sup>3</sup>)<br>FOS-DECAFE 91</th><th colspan="2">Emission des savanes 10<sup>6</sup> Tonnes</th><th colspan="3">Contribution des savanes africaines (%)</th></tr><tr><th>Afrique</th><th>Monde</th><th>Tropiques<br>feux de biomasse</th><th>Monde<br>feux de biomasse</th><th>Monde<br>toutes les sources</th></tr></thead><tbody><tr><td>Biomasse<br>brûlée</td><td></td><td>1300-2500<sup>a</sup></td><td>2500-3700<sup>a</sup></td><td></td><td></td><td></td></tr><tr><td>Carbone brûlé<br>gaz</td><td></td><td>550-1050<sup>b</sup></td><td>1050-1550<sup>b</sup></td><td>16-31</td><td>14-27</td><td>6-11</td></tr><tr><td>CO<sub>2</sub></td><td></td><td>480-920<sup>c</sup></td><td>920-1360<sup>c</sup></td><td>16-31</td><td>14-27</td><td>6-11</td></tr><tr><td>CO</td><td>61</td><td>29-56</td><td>56-83</td><td>10-19</td><td>8-16</td><td>3-5</td></tr><tr><td>CH<sub>4</sub></td><td>3,1</td><td>1,5-2,9</td><td>2,9-4,2</td><td>5-9</td><td>4-8</td><td>0,4-1</td></tr><tr><td>NMHC</td><td>7,0</td><td>3,4-6,4</td><td>6,4-9,5</td><td>16-30</td><td>14-27</td><td>3-6</td></tr><tr><td>NOy</td><td>1,91</td><td>1,1-2,1</td><td>2,1-3,0</td><td></td><td></td><td></td></tr><tr><td>NOx</td><td>1,37</td><td>0,8-1,5</td><td>1,5-2,2</td><td>11-20</td><td>9-18</td><td>2-4</td></tr><tr><td>NO<sub>2</sub></td><td>0,14</td><td>0,1-0,1</td><td>0,1-0,2</td><td></td><td></td><td></td></tr><tr><td>N<sub>2</sub>O</td><td>0,145</td><td>0,2-0,3</td><td>0,3-0,5</td><td>29-43</td><td>25-38</td><td>15-23</td></tr><tr><td>NH<sub>3</sub></td><td>0,057<sup>d</sup></td><td>0,0-0,1</td><td>0,1-0,1</td><td>0-2</td><td>0-2</td><td>0-0</td></tr><tr><td>SO<sub>2</sub></td><td>4,7 x 10<sup>-2</sup></td><td>0,06-0,12</td><td>0,12-0,17</td><td>2-5</td><td>2-4</td><td>0-0</td></tr><tr><td>COS</td><td>1,7 x 10<sup>-3</sup></td><td>0,02-0,04</td><td>0,04-0,06</td><td>25-50</td><td>22-44</td><td>1-3</td></tr><tr><td>O<sub>3</sub></td><td>30<sup>e</sup></td><td>270</td><td>410</td><td>75</td><td>64</td><td>25</td></tr><tr><td>Aérosols</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Cf</td><td>9</td><td>5,0-9,5</td><td>9,5-14,0</td><td>25-48</td><td>20-38</td><td>7-13</td></tr><tr><td>NCN</td><td>2 x 10<sup>11</sup>g</td><td>2,6 10<sup>13</sup><br>- 5,0 10<sup>13</sup></td><td>5,0 10<sup>13</sup><br>-7,4 10<sup>13</sup></td><td>16-31</td><td>14-27</td><td>-</td></tr><tr><td>IN<sup>h</sup></td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td></tr></tbody></table>

<sup>a</sup> Toutes les données sont exprimées en 10<sup>6</sup> tonnes de matière sèche.

<sup>b</sup> 42% du carbone.

<sup>c</sup> 88% du carbone émis comme CO<sub>2</sub>; émission exprimée en 10<sup>6</sup> tonnes par an C(CO<sub>2</sub>), N(NOx), S(SO<sub>2</sub>), etc.

<sup>d</sup> NH<sub>3</sub> en phase flaming.

<sup>e</sup> émissions globales et ozone d'après *Andreae*, 1991.

<sup>f</sup> aérosols carbonés < 1mm, gC/gC (matière sèche) seulement en phase flaming; émission globale d'après Cachier, 1992.

<sup>g</sup> noyaux de condensation nuageuse par gramme de matière brûlée (0,1% de saturation).

<sup>h</sup> noyaux glacogènes.

<caption>Tableau 1. Taux d'émission des gaz en traces et des aérosols; production de NCN dans les feux de savane; émissions annuelles des savanes, et contribution des savanes africaines aux émissions globales (d'après J.P. Lacaux et al., [1993]).</caption>

écosystèmes naturels ou perturbés, et la combustion de biomasse, principalement les feux de savanes, qui sont eux, avant tout d'origine anthropique.

Les processus biogéniques produisent principalement des émissions d'hydrocarbures par la végétation, de composés azotés par les sols, et de méthane dans les zones humides naturelles. Il a été montré que les feux de végétation produisent, dans des proportions variées, presque tous les gaz traces réactifs ou radiativement actifs, et des aérosols. Cependant, du fait de la composition chimique de la phytomasse, les deux types d'émissions, naturelles ou par les feux, sont surtout constitués de composés du carbone et de l'azote. En conséquence, les composés à l'état de traces émis par ces deux sources sont pratiquement les mêmes.

## 2.1. Émissions naturelles

### 2.1.1. Émissions de composés organiques.

Les plantes émettent un grand nombre de composés organiques volatils (COV), et les émissions d'isoprène, de mono et sesquiterpènes, et de composés organiques oxygénés comme les alcools, aldéhydes et cétones dépendent dans un large mesure des espèces de végétation. La canopée des forêts tropicales est la principale source d'isoprène [Gregory et al., 1986; Zimmerman et al., 1988]. Dans le cas de l'isoprène, qui n'est pas stocké par les feuilles, l'émission présente une très forte dépendance à la lumière et une dépendance à la température similaire à celle du

transport photosynthétique des électrons. La production de terpénoides dépend également, des processus photosynthétiques, mais des quantités importantes sont stockées dans la structure de certaines plantes et l'émission n'est par conséquent pas couplée à la photosynthèse et en principe non affectée par la lumière. Une partie des composés organiques émis par la végétation tropicale, en particulier l'isoprène, subit très rapidement des dégradations photo-oxydantes qui peuvent conduire à la production d'un aérosol carboné biogénique submicronique [Cachier et al., 1993] L'émission d'isoprène (420 MT[C]/an; MT = Million de Tonnes) représente environ un tiers des émissions totales de COV, elle est comparable à l'émission globale de méthane.

### 2.1.2. Émissions de méthane

Le méthane est produit par des bactéries méthanogènes strictement anaérobies. Celles-ci se développent dans des milieux anoxiques naturels tels les zones marécageuses des hautes latitudes ou des régions équatoriales, la panse des termites, mais aussi dans des milieux anoxiques générés directement ou indirectement par l'homme: rizières, panse des ruminants, décharges d'ordures ménagères ou industrielles. Au contraire la plupart des sols aérés sont des puits de méthane atmosphérique grâce à l'oxydation bactériologique du méthane par les bactéries méthanotrophes.

En milieu tropical ou les sources anthropogéniques liées à l'exploitation, au transport et à l'utilisation du gaz naturel sont faibles, les sources dominantes sont les zones inondées
<figure><img src="image_127.png" /><figcaption>Figure 1. Bilan annuel des sources de méthane en Afrique (Les parties hachurées représentent l'incertitude de l'évaluation (d'après Delmas et al., [1992])</figcaption></figure>

naturelles, la culture du riz et l'elevage de ruminants, les deux dernières s'accroissant avec la population.

Les zones inondée naturelles émettent environ 100 MT [CH₄]/an, c'est à dire 20% des émissions globales. Les forêt inondées et marécages des zone tropicales représentent environ le tiers de ces émissions. Les flux de méthane de l'ordre de 100 mg [CH₄]/m²/jour ont été observés dans le bassin de l'Amazone et du Congo [Bartlett et al., 1988 ; Tathy et al., 1992]. La culture du riz est pratiquée principalement dans les régions tropicales. Les bilans globaux incluent généralement des émissions de CH₄ par les rizières de l'ordre de 100MT/an [Cicerone et Oremland, 1988 ; Khalil et Rasmussen, 1990 ; Fung et al., 1991]. On estime que la culture du riz contribue pour 230 ppbv au stock atmosphérique actuel de CH₄ dans l'atmosphère, et pourrait être responsable d'environ 20 à 30% de l'accroissement de méthane au cours du siècle écoulé [Khalil et Rasmussen, 1991]. Les ruminants sont la troisième source majeure de méthane (80 MT/an); cependant elle n'est pas spécifique aux tropiques. La distribution géographique de cette source décrite par Lerner et Matthews [1987] est liée à la fois aux densités de population animale et aux types de régimes alimentaires [Johnson et al., 1993]. Au contraire, les émissions de méthane par les termites se localisent principalement dans les régions tropicales; cependant cette source est maintenant classée comme une source mineure de méthane, récemment estimée à 27 MT/an [Judd et al., 1993].

Plus spécifiquement en Afrique tropicale, une étude effectuée dans le cadre du programme DECAFE [Tathy et al., 1992] a montré que l'émission de méthane y est dominée par trois sources qui sont les zones inondées, les ruminants et la combustion de la biomasse. Ces trois sources, sensiblement équivalentes représentent chacune environ 10 Tg(CH₄)/an. (Figure 1). Les émissions totales pour le continent Africain sont évaluées entre 20 et 40 MT/an [Delmas et al., 1992], soit environ 1/10 de l'ensemble des sources connues; 30% du méthane produit serait d'origine abiotique puisque issu de combustions, ce qui est spécifique à ce continent (Tableau 1).

### 2.1.3. Les émissions de composés azotés par les sols.

Les oxydes d'azote (NOₓ) jouent un rôle déterminant dans la chimie de l'atmosphère, le cycle catalytique de conversion de monoxyde d'azote (NO) en dioxyde d'azote (NO₂) étant désormais la première source d'ozone troposphérique; il est donc particulièrement important de comprendre et quantifier les sources et les puits principaux. Les émissions biogéniques au niveau des sols résultent de deux processus cruciaux pour le cycle de l'azote dans les écosystèmes naturels; la nitrification (conversion de l'ammonium (NH₄⁺) en nitrate (NO₃⁻)) et la dénitrification (la réduction subséquente de NO₃⁻ en diazote (N₂)), chacun d'eux étant des processus bactériens. Le monoxyde d'azote et l'oxyde nitreux (N₂O) sont des produits intermédiaires ou dérivés de ces deux séquences de réactions. Plusieurs communautés de bactéries aérobie et anaérobie sont impliquées dans le cycle biologique complexe de l'azote dans les sols. Une production abiotique de NO résultant de la décomposition chimique des nitrites peut aussi se produire, en fonction du pH et du contenu en matière organique du sol [Blackmer et Cerrato, 1986 ; Serça et al., 1993]. Seule la dénitrification est reconnue comme un puits biologique significatif pour NO et N₂O.

Les taux de nitrification et dénitrification et le dégagement concomitant de NOₓ sont extrêmement hétérogènes dans l'espace et dans le temps, et l'élucidation des mécanismes qui contrôlent les émissions naturelles de NOₓ est un thème de recherche actif. En général, la réserve du substrat (NH₄⁺ pour la nitrification et NO₃⁻ pour la dénitrification) commande le taux de ces émissions, mais celles-ci dépendent aussi de l'environnement chimique et physique du sol (c'est à dire la texture du sol, la température, la teneur en eau, la teneur en oxygène, le pH, l'apport en carbone). Au niveau de l'écosystème, un certain nombre de régulations d'ordre plus élevé entrent en jeu, dont les taux de recyclage de l'azote, la production primaire nette, et la répartition carbone/azote dans la biomasse des plantes. Dans un écosystème non perturbé, l'azote limite presque toujours la production de plantes et il existe des mécanismes qui réduisent les pertes du système. Les perturbations de l'écosystème, qu'elles soient naturelles ou anthropiques (feu, pâturage, récolte, fertilisation) tendent à rompre ces équilibres naturels, ce qui accroît la volatilisation de l'azote du sol et son injection dans la troposphère.

Cette variabilité potentiellement élevée des émissions de composés azotés des sols, liée notamment aux perturbations anthropiques a été observée dans divers écosystèmes en Afrique [Serça et al., 1994 ; Serça, 1995], ainsi qu'en Amérique du Sud [Sanhueza et al., 1990].

### 2.2. Émission par combustion de biomasse

Outre les feux naturels, la combustion de biomasse due aux activités humaines est certainement une des plus anciennes pollutions anthropiques de l'atmosphère. La combustion de matière organique constitue un processus d'oxydation qui produit essentiellement de la vapeur d'eau et du dioxyde de carbone, dans des conditions idéales de combustion complète, d'après la réaction simplifiée :

$$[CH_2O] + O_2 \rightarrow CO_2 + H_2O$$

où [CH₂O] représente la composition moyenne de la matière organique (hydrate de carbone). Dans les feux naturels, la réserve en oxygène n'est jamais suffisante, par conséquent la combustion est incomplète et provoque la formation de composés réduits comme le méthane. Pour la plupart des constituants mineurs atmosphériques, les émissions des feux de biomasse dépendent de deux paramètres: la quantité de matière végétale brûlée, et le facteur d'émission (EF). Ce facteur est
<figure><img src="image_128.png" /><figcaption>Figure 2. Biomasses brulées annuellement dans le monde. Les zones ombrées représentent l'incertitude (valeur la plus basse et valeur la plus haute pour un type de combustion). La biomasse totale brulée est estimée en moyenne à 6,8 Gt de matière sèche (d'après Ward et Hao, [1992])</figcaption></figure>

généralement calculé en référence au CO₂ (ou à la somme CO₂+CO), sur la base du taux d'émission expérimental pour un composé donné X : EF (X) = δ[X]/δ[CO₂] (correspondant à l'excès de ces deux constituants dans les panaches de feu relativement à leur concentration moyenne dans l'air ambiant). Le facteur d'émission dépend largement du type de biomatériau brûlé et des conditions de brûlage, qui déterminent l'importance relative des deux phases principales de combustion, le flaming (combustion avec flammes) et le smoldering (combustion sans flammes proche de la pyrolyse).

### 2.2.1. Distribution des sources de combustion de biomasse dans les zones tropicales.

Parmi les différentes activités de combustion que l'on trouve dans les régions tropicales, presque toutes sont d'origine anthropique. Les principaux facteurs commandant ces activités sont liés au besoin, pour les pays en voie de développement, de sources d'énergies accessibles à des populations croissantes. Le défrichement pour l'utilisation permanente des terres pour l'agriculture et l'élevage, l'agriculture itinérante et l'utilisation de bois de chauffage sont des exemples de telles activités. Les feux de savanes, quant à eux, sont plus liés à des pratiques tradition-

nelles qu'économiques. C'est un phénomène très ancien et très répandu. Les régions de savanes brûlées représentent environ 7,5 millions de km² par an [Hao et al., 1990]. En Afrique, plus de 70% des zones de savane sont potentiellement soumises au feu chaque année. En conséquence l'Afrique est le continent qui contribue le plus à la pollution de l'atmosphère par les combustions de biomasse [Hao et al., 1990; Cachier, 1992]. Les feux se rencontrent surtout en savane humide et sont présents selon des cycles saisonniers avec un maximum d'intensité au milieu de la saison sèche, janvier-février dans l'hémisphère Nord, août-septembre dans l'hémisphère Sud [Delmas et al., 1991b; Cahoon et al., 1992].

Le défrichement de la forêt n'affecte pas la forêt tropicale humide avec la même ampleur sur les différents continents. Chaque année, on estime que 0,6% des régions forestières sont détruites et brûlées, avec actuellement la moitié de la déforestation localisée en Amérique du Sud. La déforestation dans les régions tropicales a commencé à la fin du siècle dernier mais elle a considérablement augmenté pendant les trente dernières années à la suite de l'accroissement de population. En Afrique de l'Ouest, par exemple, 70% des régions forestières existant au début du siècle ont déjà été détruites.
L'agriculture itinérante - une pratique cyclique de culture et jachère sur le même terrain - est souvent accompagnée de combustion de biomasse. On estime que l'agriculture itinérante est pratiquée par 200 millions de personnes dans le monde, principalement dans les régions tropicales, avec un défrichement annuel de 20 à 60 millions d'hectares [Seiler et Crutzen, 1980; Hao et al., 1990].

La collecte de bois de chauffe affecte les écosystèmes forestiers comme les savanes. Le bois de chauffe est la principale source d'énergie pour l'utilisation domestique dans la plupart des pays tropicaux et cela soulèvera un problème économique important dans les dix prochaines années : l'augmentation de la demande en bois de chauffe entraîne une diminution des ressources en bois dans de nombreux pays. Le bois de chauffe peut être utilisé directement pour la cuisine ou le chauffage, ou peut être transformé en charbon de bois principalement pour la fourniture en énergie domestique des grandes cités. La production de charbon augmente ainsi proportionnellement à l'accroissement urbain dans les pays en voie de développement.

La dernière source de combustion de biomasse dans les tropiques est le brûlage saisonnier des résidus agricoles. C'est le cas pour la paille de riz en Asie du Sud-Est ou la canne à sucre qui est presque entièrement brûlée sur le champ avant la récolte et peut être une source locale très importante de pollution [Kirchhoff et al., 1989].

### 2.2.2 La biomasse brûlée annuellement à travers le monde

Une récente estimation de la biomasse brûlée annuellement dans le monde est présentée sur la figure 1. La quantité de biomasse brûlée dans un écosystème est généralement calculée d'après l'équation suivante :

$$M = A \times B \times a \times b \quad [Seiler et Crutzen, 1980]$$

où A est la surface soumise au feu, B est la densité de biomasse totale, a est la proportion de biomasse épigée par rapport à la biomasse totale, et b est le facteur de combustion qui représente la fraction de biomasse épigée réellement brûlée.

Un certain nombre d'investigations, basées sur des observations, ont été conduites dans différents écosystèmes pour estimer ces paramètres. Elles concernent la biomasse réellement soumise au feu et les différentes caractéristiques des combustions [Hao et al., 1990; Menaut et al., 1991; Kauffman et al., 1992; Ward et al., 1992].

### 2.2.3 Emissions de gaz traces et d'aérosols par combustion de biomasse

Les émissions de gaz traces et d'aérosols par les feux de biomasse dépendent de la quantité de matériau brûlé et du facteur d'émission (EF) de chaque composé. La quantification de l'émission doit prendre en compte les caractéristiques physiques de la végétation brûlée, qui déterminent l'importance relative des deux phases principales de la combustion (flaming et smoldering). En conséquence la proportion des différents gaz traces ou aérosols émis pendant la combustion dépend fortement du type de combustion et des matériaux brûlés. On a montré par exemple que le facteur d'émission du méthane ou des aérosols peut varier d'un facteur supérieur à 10 entre feux de savane, feux de forêt, et production de charbon de bois [Delmas et al., 1991a; Cachier et al., 1994]. Pour les composés de l'azote, les facteurs d'émissions dépendent aussi largement de la composition chimique des plantes car la teneur en azote des plantes peut varier selon le type de plante, et aussi selon la saison pour un type de plante donné (Figure 2).

Les aérosols produits par les feux de biomasse contiennent principalement du carbone sous forme d'espèces organiques

polymérisées parmi lesquelles la matière la plus réfractaire est constituée par le carbone-suie (carbone élémentaire). La teneur des particules en carbone suie (3 à 25%) est fortement variable selon, en premier lieu, les paramètres physiques de la combustion [Cachier et al., 1991]. Les particules contiennent également tout un cortège d'éléments traces qui proviennent de la matrice végétale mais aussi, pour une part importante, de la remobilisation de poussières d'origines diverses déposées sur les végétaux. Il semble cependant que la richesse en potassium des aérosols (rapport potassium/carbone-suie~0,4) soit une caractéristique des émissions des feux de végétation et permettent de les distinguer des aérosols industriels. Les particules atmosphériques affectent le bilan radiatif de la Terre, à la fois directement à cause du rayonnement solaire qu'elles rétrodiffusent et indirectement en se comportant comme noyaux de condensation (NCN) qui modifient la couverture nuageuse. Des estimations sur la production de gaz traces, de particules carbonées et de NCN par les feux de savane en Afrique sont données dans le tableau 1 [d'après Lacaux et al., 1993b] qui illustre l'importance régionale et globale de ces sources.

## 3. INFLUENCE DES ÉMISSIONS DE CONSTITUANTS MINEURS EN ZONE TROPICALE SUR LA CHIMIE TROPOSPHÉRIQUE

### 3.1. Processus atmosphériques

Les composés en traces émis par la biosphère au travers de processus naturels ou combustifs sont à la base des phénomènes photochimiques, qui en présence de lumière solaire déterminent le cycle des oxydants dans la troposphère. Les processus de chimie atmosphérique dans les régions tropicales sont gouvernés par plusieurs facteurs :

(i) La météorologie qui détermine les régimes de pluviosité, la circulation régionale des masses d'air, et les transferts verticaux de la surface à la couche limite et l'atmosphère libre; tous ces processus régulent l'interaction entre les deux principaux écosystèmes : la forêt équatoriale et la savane.

(ii) La présence de sources et puits de composés en trace détaillés précédemment.

(iii) La chimie et la photochimie dans l'atmosphère qui mènent, pour les composés réactifs, à la transformation de produits primaires instables en produits secondaires plus stables qui peuvent être transportés dans les masses d'air loin des régions sources et donc affecter la composition chimique de l'atmosphère à l'échelle régionale et globale. Ces transformations s'effectuent au travers de processus homogènes et/ou hétérogènes. Les processus hétérogènes sont particulièrement marqués dans la zone équatoriale à cause de l'importance de la couverture nuageuse et de l'activité convective. La convection nuageuse facilite aussi les échanges verticaux de la couche limite vers la haute troposphère et donc le transport à longue distance de composés en traces issus de la biosphère tropicale continentale.

Ces processus sont illustrés sur la figure 3 représentant un schéma des interactions biosphère-atmosphère en Afrique tropicale pendant la saison sèche. Les panaches de fumée sont entraînés par l'alizé continental sec du Nord-Est (Harmattan) vers la forêt équatoriale. Les panaches individuels deviennent diffus en quelques dizaines de kilomètres, formant un voile de fumée uniformément mélangé. Les émanations donnent des oxydes d'azote, du monoxyde de carbone, des hydrocarbures réactifs qui produisent de l'ozone, ainsi que des aérosols. Les composés émis subissent des réactions atmosphériques initiales de panaches concentrés, suivies par un mélange dans une grande partie de la troposphère et un retour possible dans la couche limite. La masse d'air, contaminée par les produits de combustion des feux de
<figure><img src="image_129.png" /><figcaption>Figure 3. Emission d'oxydes d'azote par les feux de savane: relation (linéaire) entre le facteur d'émission et la teneur en azote des plantes.</figcaption></figure>

savanes, est piégée entre la couche de mousson et le flux d'Ouest permanent. Au contraire, pendant la saison humide, la zone de savane est influencée par la mousson humide du Sud-Est qui traverse la forêt tropicale, laquelle est à son tour influencée par l'air océanique. Les composés azotés émis par la forêt et transportés vers la savane constituent des nutriments essentiels pour la croissance végétale. Ceci illustre la complexité des échanges biosphère-atmosphère, et la chimie troposphérique qui leur est associée dans les régions tropicales.

## 3.2. Cycle des oxydants : formation de l'ozone

Ces dix dernières années, les instruments satellitaires, les programmes d'observations intensives et les réseaux-sol de surveillance, qui ont fourni des données sur l'ozone dans la ceinture tropicale, ont mis en évidence l'impact des feux sur la chimie atmosphérique. Les données-satellite ont été obtenues à partir des instruments TOMS (Total Ozone Mapping Spectrometer) et SAGE (Stratospheric Aerosol Gas Experiment) entre 1979 et 1989. Les taux d'ozone troposphérique ont été calculés en soustrayant le taux d'ozone dans la stratosphère, tiré des mesures de SAGE, de la mesure concomitante d'ozone total acquise par TOMS au même endroit. La différence est d'environ 10% et correspond à l'ozone troposphérique [Fishman et al., 1990]. Les données troposphériques ont montré un maximum bien défini d'ozone supérieur à 45 unités Dobson [1 D.U. = 2.69 x 10¹⁶ molécules d'ozone par cm²] au large de la côte Ouest de l'Afrique en septembre et octobre (Figure 4). Cette période correspond au maximum de densité de feux observés dans la zone de savane au Sud de l'Afrique. Les sondages d'ozone effectués entre 1990 et 1992, dans trois stations - à Brazzaville (Congo), a l'île de l'Ascension et à Natal (Brésil) - dans le cadre de l'expérience PréTRACE ont confirmé l'hypothèse selon laquelle l'Ozone produit au dessus de l'Afrique du Sud est transporté par les flux d'Est des basses couches vers l'océan Atlantique [Cros et al., 1992]. Une partie de l'ozone produit dans la couche limite est introduit dans la moyenne et la haute troposphère où les vents dans les régions tropicales et subtropicales sont généralement d'Ouest, expliquant l'augmentation des concentrations d'ozone observée dans l'océan Indien entre l'Afrique et l'Australie [Fishman et al., 1991].

Les observations aéroportées dans les régions tropicales aboutissent à des conclusions similaires. Dans les régions affectées par les feux de biomasse tropicale, les concentrations d'ozone se situent entre 50 et 100 ppb. En outre les profils verticaux présentent une structure typique en couches discrètes à forte concentration d'ozone entre 1 et 5 km d'altitude (Figure 5) [Marenco et al., 1990]. Cette structure verticale en couches indique l'importance des transports convectifs pour la distribution des produits primaires et secondaires des feux de biomasse.

En plus des émissions associées aux feux de biomasse, les émanations d'espèces volatiles par les sols et la végétation contribuent directement à la production photochimique d'espèces chimiquement réactives. A proximité des zones de transition entre savane et forêt, des teneurs en ozone atteignant 70 ppbv sont couramment observés entre 1 et 3 km d'altitude (Figure 6). Les émissions de feux de biomasse et les produits secondaires de la photochimie associée atteignent les régions de forêt tropicale [Andreae et al., 1992b] et peuvent réagir avec les composés en traces émis par la forêt, ce qui provoque l'augmentation des niveaux d'ozone dans la basse atmosphère [Lopez et al., 1992].

## 3.3. Impact radiatif

On sait que les feux de biomasse sont une des principales causes de l'effet de serre additionnel et de sa conséquence possible: le réchauffement global de l'atmosphère. Quel que soit le type de combustion, le principal composé émis, le gaz carbonique (CO₂), représente entre 85 et 90 % du carbone injecté dans l'atmosphère [Lobert et al., 1991]. Cependant, comme l'indiquent plusieurs articles [Crutzen et Andreae, 1990 ; Andreae, 1991], cette émission ne correspond, pas toujours, à un flux net de carbone vers l'atmosphère. Les feux de savanes, de même que ceux des friches agricoles ne contribuent pas à une augmentation à long terme du dioxyde de carbone dans l'atmosphère. En effet le CO₂ émis par de telles combustions est à nouveau ôté du réservoir atmosphérique et réincorporé dans la biomasse par l'activité photosynthétique pendant la repousse de la végétation dans les semaines ou les mois qui suivent les feux. Seuls les feux de forêt pour le défrichage de nouvelles terres et l'utilisation de bois comme combustible prennent directement
<figure><img src="image_130.png" /><figcaption>Figure 4. Schéma conceptuel des échanges biosphère-atmosphère et des processus photochimiques associés, à l'interface forêt-savane, en saison sèche, en Afrique tropicale.</figcaption></figure>

<figure><img src="image_131.png" /><figcaption>Figure 5. Répartition de l'ozone dans la colonne troposphérique (en unité Dobson) entre 50°N et 50°S obtenue par détermination satellitaire du résidu troposphérique</figcaption></figure>
<figure><img src="image_132.png" /><figcaption>Figure 6. Profils verticaux de O₃, CO, CO₂, NOy mesurés pendant l'expérience Tropoz I au dessus de la Côte d'Ivoire, en saison sèche (d'après Marenco et al. [1990]).</figcaption></figure>

part dans l'effet de serre additionnel. Néanmoins, tous les processus de combustion de biomasse constituent une source nette pour presque tous les gaz en traces autres que le CO₂, de même que pour les aérosols. La dégradation bactériologique de matière organique sous des conditions aérobies aboutit essentiellement à l'émission de CO₂ alors que la combustion de biomasse est accompagnée par l'émission ou la formation secondaire d'autres composés en traces radiativement actifs (CH₄, O₃). D'autre part, les feux de biomasse produisent également des aérosols qui pourraient avoir des effets inversent de ceux des gaz à effet de serre, à cause de leur propriétés optiques et de leur rôle de noyaux de condensation nuageuse. De plus la déforestation peut aussi aboutir à un destockage relativement rapide du carbone des sols. Globalement la déforestation dans les tropiques est supposée contribuer pour plus de 20% au réchauffement anthropique global [IPCC, 1990].

Les forêts tropicales constituent également la principale source de protoxyde d'azote atmosphérique, estimée à environ 2MT d'azote par an. Le protoxyde d'azote est un gaz à effet de serre et il intervient également dans la chimie de l'ozone stratophérique. Il a été démontré que la conversion des forêts en systèmes agricoles tend à perturber le recyclage des éléments nutritifs et donc à modifier les émissions de protoxyde d'azote, au moins pendant quelques années après la déforestation [Keller et al., 1993].

### 3.4. Acidification des pluies

Des taux élevés de précipitations acides sont signalées dans des régions tropicales éloignées [Lacaux et al., 1992b]. Pendant

<figure><img src="image_133.png" /><figcaption>Figure 7. Comparaison de la distribution verticale de l'ozone observée pendant l'expérience DECAFE 88 avec les profils moyens du Pacifique équatorial et du bassin amazonien central pendant les saisons sèches et humides. La zone ombrée définit une enveloppe comprenant toutes les données de l'expédition DECAFE 88 (saison sèche au nord du Congo) (d'après Andreae et al. [1992])</figcaption></figure>