Frédéric Frappart
email: frederic.frappart@legos.obs-mip.fr

Pierre Zeiger
email: pierre.zeiger@legos.obs-mip.fr

Julie Betbeder
email: julie.betbeder@cirad.fr

Valéry Gond
email: valery.gond@cirad.fr

Régis Bellot
email: regis.bellot@orange.fr

Nicolas Baghdadi
email: nicolas.baghdadi@teledetection.fr

Fabien Blarel
email: fabien.blarel@legos.obs-mip.fr

José Darrozes
email: jose.darrozest@get.omp.eu

Luc Bourrel
email: luc.bourrel@ird.fr

Frédérique Seyler
email: frederique.seylert@ird.fr

remote sensing

Keywords: radar altimetry, wetlands, surface water storage, Congo

Surface water storage in floodplains and wetlands is poorly known from regional to global scales, in spite of its importance in the hydrological and the carbon balances, as the wet areas are an important water compartment which delays water transfer, modifies the sediment transport through sedimentation and erosion processes, and are a source for greenhouse gases. Remote sensing is a powerful tool for monitoring temporal variations in both the extent, level, and volume, of water using the synergy between satellite images and radar altimetry. Estimating water levels over flooded area using radar altimetry observation is difficult. In this study, an unsupervised classification approach is applied on the radar altimetry backscattering coefficients to discriminate between flooded and non-flooded areas in the Cuvette Centrale of Congo. Good detection of water (open water, permanent and seasonal inundation) is above 0.9 using radar altimetry backscattering from ENVISAT and Jason-2. Based on these results, the time series of water levels were automatically produced. They exhibit temporal variations in good agreement with the hydrological regime of the Cuvette Centrale. Comparisons against a manually generated time series of water levels from the same missions at the same locations show a very good agreement between the two processes (i.e., RMSE ≤ 0.25 m in more than 80%/90% of the cases and R ≥ 0.95 in more than 95%/75% of the cases for ENVISAT and Jason-2, respectively). The use of the time series of water levels over rivers and wetlands improves the spatial pattern of the annual amplitude of water storage in the Cuvette Centrale. It also leads to a decrease by a factor of four for the surface water estimates in this area, compared with a case where only time series over rivers are considered.

Introduction

Floodplains and wetlands cover at least 12.1 × 10 6 km 2 (~8%) of the land surfaces of the Earth [START_REF] Davidson | Global extent and distribution of wetlands: Trends and issues[END_REF][START_REF] Mitsch | Wetland Ecosystems[END_REF]. They play a major role in the water cycle through river flow variability, flood mitigation, groundwater recharge and water quality improvement [START_REF] Mitsch | Wetland Ecosystems[END_REF][START_REF] Bullock | The role of wetlands in the hydrological cycle[END_REF][START_REF] Acreman | How wetlands affect floods[END_REF]. They were identified as one of the most productive ecosystems as well as a major contributor of biodiversity within a landscape [START_REF] Ghermandi | The Economic Value of Wetland Conservation and Creation: A Meta-Analysis[END_REF][START_REF] Junk | The comparative biodiversity of seven globally important wetlands: A synthesis[END_REF][START_REF] Mitra | An appraisal of global wetland area and its organic carbon stock[END_REF][START_REF] Webb | Effects of Local and Landscape Variables on Wetland Bird Habitat Use During Migration Through the Rainwater Basin[END_REF][START_REF] Haukos | The importance of playa wetlands to biodiversity of the Southern High Plains[END_REF][START_REF] Gibbons | Terrestrial habitat: A vital component for herpetofauna of isolated wetlands[END_REF]. They also have an important role in the global carbon cycle as 16 to 33% of the soil carbon pool is stored in the wetlands [START_REF] Maltby | Carbon dynamics in peatlands and other wetland soils regional and global perspectives[END_REF][START_REF] Mitsch | Wetlands, carbon, and climate change[END_REF] and as 20 to 25% of the methane emissions originated from the wetlands [START_REF] Whalen | Biogeochemistry of Methane Exchange between Natural Wetlands and the Atmosphere[END_REF][START_REF] Bergamaschi | Satellite charactography of atmospheric methane from SCIAMACHY on board ENVISAT: 2. Evaluation based on inverse model simulations[END_REF][START_REF] Bloom | Large-scale controls of methanogenesis inferred from methane and gravity spaceborne data[END_REF][START_REF] Ringeval | An attempt to quantify the impact of changes in wetland extent on methane emissions on the seasonal and interannual time scales[END_REF].

In spite of their importance, the water stored in floodplains and wetlands and its temporal variations are poorly known from regional to global scales. Before the launch of NASA/CNES Surface Water and Ocean Topography (SWOT) in 2022, which will provide surface water elevation over inland water bodies at a spatial resolution of 100 m [START_REF] Biancamaria | The SWOT Mission and Its Capabilities for Land Hydrology[END_REF], anomalies of surface water storage are currently derived by: (i) estimating water level changes using an interferometric synthetic aperture radar (SAR) (InSAR) [START_REF] Alsdorf | Amazon floodplain water level changes measured with interferometric SIR-C radar[END_REF][START_REF] Lee | Mapping wetland water depths over the central Congo Basin using PALSAR ScanSAR, Envisat altimetry, and MODIS VCF data[END_REF], (ii) filling a digital elevation model (DEM) with surface water extent estimates from remotely sensed observations through a hypsometric curve [START_REF] Papa | Surface freshwater storage and variability in the Amazon basin from multi-satellite observations, 1993-2007[END_REF][START_REF] Salameh | Fifteen years (1993-2007) of surface freshwater storage variability in the ganges-brahmaputra river basin using multi-satellite observations[END_REF], (iii) combining surface extent products derived from satellite images with radar altimetry based water levels [START_REF] Frappart | Floodplain water storage in the Negro River basin estimated from microwave remote sensing of inundation area and water levels[END_REF][START_REF] Frappart | Surface freshwater storage and dynamics in the Amazon basin during the 2005 exceptional drought[END_REF], or (iv) solving the water balance equation combining various remotely sensed observations (i.e., anomalies of terrestrial water storage (TWS) from the Gravity Recovery And Climate Experiment (GRACE), water levels from radar altimetry, rainfall from Global Precipitation Climatology Project (GPCP), shuttle radar terrain model (RTM), digital elevation model (DEM), synthetic aperture radar (SAR) from Japanese Earth Resources -1 (JERS-1) and multi-spectral from MODerate resolution imaging spectroradiometer (MODIS) images) [START_REF] Alsdorf | Seasonal water storage on the Amazon floodplain measured from satellites[END_REF][START_REF] Lee | Characterization of terrestrial water dynamics in the Congo Basin using GRACE and satellite radar altimetry[END_REF]. The latter approach needs the creation of a network of altimetry-based stations of water levels or virtual stations (VSs). Currently, these VSs are mostly obtained through a manual process using dedicated softwares such as the Virtual ALtimetric Stations (VALS) software [START_REF] Da Silva | Water levels in the Amazon basin derived from the ERS 2 and ENVISAT radar altimetry missions[END_REF], the Multi-mission Altimetry Processing Software (MAPS) [START_REF] Frappart | Preliminary Assessment of SARAL/AltiKa Observations over the Ganges-Brahmaputra and Irrawaddy Rivers[END_REF][START_REF] Normandin | Evolution of the performances of radar altimetry missions from ERS-2 to Sentinel-3A over the Inner Niger Delta[END_REF], or the Altimetry Time Series (AlTiS) software [START_REF] Frappart | Evaluation of the Performances of Radar and Lidar Altimetry Missions for Water Level Retrievals in Mountainous Environment: The Case of the Swiss Lakes[END_REF], except for large rivers and lakes. This processing is time consuming and does not allow a complete coverage of extensive areas such as the wetlands present in the large river basins, especially when considering a large number of satellite missions. Automatic processing to produce RA-based water levels relies on the identification of the cross-section between a permanent waterbody and the RA ground-tracks (e.g., [START_REF] Schwatke | DAHITI-An innovative approach for estimating water level time series over inland waters using multi-mission satellite altimetry[END_REF][START_REF] Okeowo | Automated Generation of Lakes and Reservoirs Water Elevation Changes from Satellite Radar Altimetry[END_REF]). Therefore, these processes are applied to cross-sections located on rivers, lakes and reservoirs.

Yet, several studies already shown that radar altimetry (RA) can be used for the monitoring of land surface properties. Spatio-temporal changes in radar altimetry backscattering coefficients (σ 0 ) were related to land cover types (ice, arid and semi-arid areas, wetlands, forests, . . . ) and hydrological changes (floods, soil moisture, . . . ) from regional to global scales [START_REF] Papa | Use of the Topex-Poseidon dual-frequency radar altimeter over land surfaces[END_REF][START_REF] Legrésy | ENVISAT radar altimeter measurements over continental surfaces and ice caps using the ICE-2 retracking algorithm[END_REF][START_REF] Papa | Inundated wetland dynamics over boreal regions from remote sensing: The use of Topex-Poseidon dual-frequency radar altimeter observations[END_REF][START_REF] Fatras | Spaceborne altimetry and scatterometry backscattering signatures at C-and Ku-bands over West Africa[END_REF][START_REF] Frappart | Radar altimetry backscattering signatures at Ka, Ku, C, and S bands over West Africa[END_REF][START_REF] Frappart | Backscattering signatures at Ka, Ku, C and S bands from low resolution radar altimetry over land[END_REF]. Over West Africa, signature of floods under the dense vegetation canopy of the Congo River Basin was observed at S, C, Ku and Ka bands through an increase in the backscattering coefficient of various altimetry missions (i.e., ENVISAT, Jason-2 and SARAL) [START_REF] Frappart | Radar altimetry backscattering signatures at Ka, Ku, C, and S bands over West Africa[END_REF]. These results are very encouraging for using the backscattering coefficient to automatically identify radar altimetry measurements over water under dense vegetation cover.

The objectives of this study are to (i) demonstrate the possibility to identify water under the altimeter ground-tracks using the backscattering coefficient, important but poorly used information derived from the RA echo [START_REF] Abdalla | Altimetry for the future: Building on 25 years of progress[END_REF], (ii) automatically create time series of water levels over inland waterbodies, including floodplains and wetlands which are not monitored, (iii) to measure the impact of the densification of the networks of altimetry-based water levels on the surface water estimates currently performed combining nadir RA and satellite images. This assessment was achieved on the Cuvette Centrale of Congo, applying an unsupervised classification technique (k-means) to radar altimetry backscattering coefficients acquired at Ku-band (ENVISAT and Jason-2) to automatically discriminate acquisitions made over water bodies (e.g., open water such as rivers, ponds, oxbow lakes, and lakes and water under vegetation as forested floodplains and wetlands) in the extensive floodplain area of the Cuvette Centrale in the Congo River basin. The optimal selection of the cluster numbers was performed using the Calinski-Harabasz criterion [START_REF] Caliñski | A Dendrite Method Foe Cluster Analysis[END_REF]. The results of the clustering technique were validated against an unsupervised classification of the Cuvette Centrale validated using several data sources [START_REF] Betbeder | Mapping of central africa forested wetlands using remote sensing[END_REF]. Based on the results of the clustering technique, the time series of water levels were generated at the cross-sections with the river network and over the floodplains to construct a network of VSs. Over rivers, the resulting time series of water level were compared, due to the lack of in situ data in this area, this study used the altimetry-based water stages from the Hydroweb database [41]. Then, water level maps were obtained combining surface water extent from PALSAR-1 images (see [START_REF] Betbeder | Mapping of central africa forested wetlands using remote sensing[END_REF] for details about the processing of the SAR images) and altimetry-based water levels over rivers and floodplains.

Study Area

The Congo River Basin (Figure 1a) is the largest drainage basin in Africa. Its shape is almost circular, covering an area of ~4 million km 2 [START_REF] Devroey | Annuaire Hydrologique du Congo et Ruanda-Urundi[END_REF]. It is the second river in the world, after the Amazon River, in terms of discharge (~40,000 m 3 •s -1 ) [START_REF] Bricquet | Les Ecoulements Du Congo a Brazzaville Et La Spatialisation Des Apports[END_REF]. Its main waterway is the Congo River, with a length of ~6650 km from its source in the southeast of the Democratic Republic of Congo to the Atlantic Ocean [START_REF] Runge | The Congo River, Central Africa[END_REF]. The Congo River Basin is mainly located in the equatorial and tropical savanna climate area according to the Köppen-Geiger climate classification [START_REF] Peel | Updated world map of the Köppen-Geiger climate classification[END_REF], but also in the humid subtropical climate in its southern part, and, marginally, in the semi-arid and desert climates in the east (Figure 1a).

The Cuvette Centrale is a vast floodplain located in the center of the Congo River Basin which extends from 3 • S to 3 • N in latitude and from 16 • E to 22 • E in longitude (Figure 1b). It covers an area of 1,176,000 km 2 with a wetland extent of 360,000 km 2 (32% of its area) [START_REF] Bwangoy | Wetland mapping in the Congo Basin using optical and radar remotely sensed data and derived topographical indices[END_REF]. It is the remnant of a lake which occupied the area during the Tertiary geological period that is now surrounded by mountains and plateaus [START_REF] Robert | Le Congo Physique[END_REF]. As this area is difficult to access, the vegetation has been poorly studied [START_REF] Campbell | The Congo river basin[END_REF]. It is mostly covered with tropical evergreen forests, and vegetation adapted to soils saturated with water (e.g., flooded forests and inundated grasslands) [START_REF] Betbeder | Mapping of central africa forested wetlands using remote sensing[END_REF][START_REF] Bwangoy | Wetland mapping in the Congo Basin using optical and radar remotely sensed data and derived topographical indices[END_REF][START_REF] Vancutsem | Mapping and characterizing the vegetation types of the Democratic Republic of Congo using SPOT VEGETATION time series[END_REF][START_REF] Verhegghen | Mapping Congo Basin vegetation types from 300 m and 1 km multi-sensor time series for carbon stocks and forest areas estimation[END_REF]. In the Cuvette Centrale are located the confluences between the Congo River and two of its major northern tributaries, the Ubangi River, and downstream, the Sangha River. The annual rainfall in the Cuvette Centrale ranges from 1400 to 1800 mm•year -1 and the potential evapotranspiration reaches 1280 mm•year -1 [START_REF]FAO World-Wide Agroclimatic Data of FAO (FAOCLIM)|Land & Water|Food and Agriculture Organization of the United Nations|Land & Water|Food and Agriculture Organization of the United Nations[END_REF]. The topography is flat with an average slope lower than 7 cm•km -1 between Kisangani and Kinshasa [START_REF] Devroey | Annuaire Hydrologique du Congo et Ruanda-Urundi[END_REF], on the 50 km upstream is the outlet of the Kasai River [START_REF] Devroey | [END_REF] and on the 450 km upstream, is the downstream part of the Ubangi River [START_REF] Yayer | Caractéristiques Hydrographiques de l'Oubangui; Comité Hydrographique du Bassin Congolais[END_REF]. The Cuvette Centrale is covered with sandy lacustrine Quaternary sediments where are located the forested wetlands [START_REF] Moukolo | Transport en solution et en suspension par le fleuve congo (zaïre) et ses principaux affluents de la rive droite[END_REF][START_REF] Censier | Caractérisation de processus d'érosion régressive par analyse sédimentologique comparée des sables du chenal et des barres du cours inférieur de l'Oubangui (République Centrafricaine, Congo, Zaïre)[END_REF]. These forested wetlands are subject to the flood pulse of the Congo River [START_REF] Junk | The flood pulse concept in river-floodplain-systems[END_REF]. Water levels of the Congo River are characterized by a bimodal flooding pattern, with a main peak of high water in November-December, and a secondary one in April-May, a more pronounced low water stage in August and one of lower intensity occurring in February-March in the Cuvette Centrale. The hydrological regime of the Ubangi (main tributary in the north of the Cuvette Centrale) is characterized by a high water period peaking in November and low water levels in March [START_REF] Wesselink | Randriamiarisoa Les régimes hydroclimatiques et hydrologiques d'un bassin versant de type tropical humide: L'Oubangui (République Centrafricaine)[END_REF]. Permanently inundated areas are located alongside the rivers. Other areas are inundated either during both flood events (discontinuously during a long time period of 6-7 months) or only during the maximum of the largest inundation event occurring in November-December (during a shorter time period of 2-3 months). More details about the Cuvette Centrale, its vegetation and hydrological regime can be found in [START_REF] Betbeder | Mapping of central africa forested wetlands using remote sensing[END_REF].

Datasets

Radar Altimetry Data

Data acquired by two altimetry missions (ENVISAT and Jason-2) were considered in this study. These two missions were operating on their nominal orbits from March 2002 to October 2010 and from June 2008 to September 2016 respectively. These two missions were chosen because they were operating when PALSAR-1 onboard ALOS-1 was in operation (see Section 3.2). A complete description of these two altimetry missions can be found in [START_REF] Benveniste | The radar altimetry mission: RA-2, MWR, DORIS and LRR[END_REF] and in [START_REF] Lambin | The OSTM/Jason-2 Mission[END_REF]. In this study, the following parameters, available at high frequency (i.e., at 18 Hz for ENVISAT and 20 Hz for Jason-2, that is to say ~350 m along the altimetry ground-tracks), were used: the time and geographical locations of acquisition, the distance between the satellite and the surface or range, the orbit of the satellite, the corrections applied to the range to account for the delays caused by the path through the atmosphere (ionosphere, dry and wet atmosphere corrections) and some geophysical effects (solid earth and pole tides), the geoid model and the backscattering coefficient (σ 0 ) at Ku-band. Radar altimeter ranges and backscattering coefficients at Ku-band are derived from the offset center of gravity (OCOG) retracking algorithm [START_REF] Wingham | New Techniques in Satellite Altimeter Tracking Systems[END_REF], found to be well-suited for land surface studies in several previous studies (e.g., [START_REF] Da Silva | Water levels in the Amazon basin derived from the ERS 2 and ENVISAT radar altimetry missions[END_REF][START_REF] Frappart | Preliminary results of ENVISAT RA-2-derived water levels validation over the Amazon basin[END_REF]). All these parameters, contained in the geophysical data records (GDR) of each mission, are made available by Centre de Topographie des Océans et de l'Hydrosphère (CTOH) [62], as well as indexes of the normalized tracks [START_REF] Blarel | Altimetry backscattering signatures at Ku and S bands over land and ice sheets[END_REF][START_REF] Frappart | An ERS-2 altimetry reprocessing compatible with ENVISAT for long-term land and ice sheets studies[END_REF] that were used to compute statistics along the tracks (see Section 4.2). The location of the ground-tracks in the Cuvette Centrale is presented in Figure 1b.

Land Cover Map of the Cuvette Centrale

It results from the merging of two unsupervised classifications obtained using the k-means clustering technique: (i) From the enhanced vegetation index (EVI) [START_REF] Huete | Development of vegetation and soil indices for MODIS-EOS[END_REF] Four classes of forest types and four classes of hydrological status were obtained from the Thorndike index [START_REF] Thorndike | Who belongs in the family?[END_REF]. Through cross-comparisons with several external datasets, the topographic map of northern Congo, ICESat lidar data of elevation and canopy height and aerial photographs, five classes were identified attributed to: More details about the classification results can be found in [START_REF] Betbeder | Mapping of central africa forested wetlands using remote sensing[END_REF].

Altimetry-Based Time Series of Water Levels

As no in situ data of water stage are available in the study area during the altimetry period, time series of water levels derived from ENVISAT and Jason-2 data over rivers in the Cuvette Centrale were used to validate the automatic estimates of river heights derived from the unsupervised classification applied to the radar altimetry backscattering coefficients. They were obtained from the Hydroweb database [41]. In this study, 24 and 32 time series from ENVISAT (35-day repeat period) and Jason-2 (10-day repeat period), respectively, were used. The length of the cross-sections between the river and the altimetry tracks range from 100 to 150 m in the upstream parts of the Cuvette Centrale to 10 km in its downstream part.

Methods

Radar Altimetry Data Pre-Processing

Radar altimetry data were preprocessed in order to:

(i) Derive the altimeter height (h) from the parameters contained in the GDR using the following equation (e.g., [START_REF] Crétaux | Hydrological applications of satellite altimetry: Rivers, lakes, man-made reservoirs, inundated areas[END_REF]):

h = H -R -Σ(∆R propagation + ∆R geophysical ) + N, ( 1
)
where H is the height of the center of mass of the satellite above the ellipsoid estimated using the Precise Orbit Determination (POD) technique, R is the altimeter range (i.e., the distance from the center of mass of the satellite to the surface taking into account instrumental corrections equals to c∆t/2 where c is the velocity of light in the vacuum and ∆t is the two-way travel time of the electromagnetic wave emitted by the radar), ∆R propagation and ∆R geophysical are the sum of the geophysical and environmental corrections applied to the range, respectively, and N is the geoid. The corrections to apply to the range for the propagation are the following:

∆R propagation = ∆R ion + ∆R dry + ∆R wet , (2)
where ∆R ion is the atmospheric refraction range delay due to the free electron content associated with the dielectric properties of the ionosphere, ∆R dry is the atmospheric refraction range delay due to the dry gas component of the troposphere, ∆R wet is the atmospheric refraction range delay due to the water vapor and the cloud liquid water content of the troposphere. And for the geophysical effects:

∆R geophysical = ∆R solid•Earth + ∆R pole , (3)
where ∆R solid•Earth and ∆R pole are the corrections accounting for crustal vertical motions due to the solid Earth and pole tides, respectively.

(ii) Obtain along-track time series of backscattering coefficient presenting no missing data in entry of the clustering method. To do so, monthly climatologies (i.e., the average of all data from each month averaged over the whole study period) of backscattering at Ku-band for every normalized index along the altimetry groundtracks were computed.

Radar Altimetry Data Clustering

K-means clustering is a statistical technique designed to assign objects to a fixed number of groups or clusters based on the analysis of a set of specified variables so that the within-cluster sum of squares is minimum [START_REF] Macqueen | Some methods for classification and analysis of multivariate observations[END_REF]. The initial cluster centers K are selected and then iteratively refined: (i) Each object is assigned to its closest cluster in terms of distance to the center of the cluster; (ii) At every iteration, each cluster center or centroid is updated to be the average of the member of the cluster.

This process is repeated until convergence. K-means++ algorithm, which improves the quality of the final solution in terms of intra-cluster distances [START_REF] Arthur | k-means++: The advantages of careful seeding[END_REF] was used. It benefits from repeated new initial cluster centroid positions. The number of times to repeat clustering using the new initial cluster centroid positions was set to 5 and the maximum number of iterations before convergence to 100. The distance used for the k-means approach is the Euclidean distance (Ed).

In this study, the k-means approach is used to perform an unsupervised classification of the radar altimetry data using the seasonality of the backscattering at Ku-band. Data used as an input of the k-means clustering are climatological monthly variations in the backscattering at Ku-band for every normalized index along the altimetry ground-tracks [START_REF] Frappart | Backscattering signatures at Ka, Ku, C and S bands from low resolution radar altimetry over land[END_REF]. The average ( ) and standard deviation (std) of the radar altimetry backscattering coefficients (σ 0 ) are computed on every along-track grid-point as follows [START_REF] Frappart | Radar altimetry backscattering signatures at Ka, Ku, C, and S bands over West Africa[END_REF][START_REF] Frappart | Backscattering signatures at Ka, Ku, C and S bands from low resolution radar altimetry over land[END_REF]: σ 0 (dB) = 10 log 10 10 σ 0 /10 , ( 4)

std(σ 0 ) (dB) = 10 log 10   1 + std 10 σ 0 /10 10 σ 0 /10   , (5)
The optimal number of clusters was determined using the Calinski-Harabasz criterion (I CH ) and the silhouette coefficient (6 and 7, respectively). I CH or variance ratio criterion was defined as the between-cluster variance and the overall within-cluster variance [START_REF] Caliñski | A Dendrite Method Foe Cluster Analysis[END_REF]:

I CH = N S -C C -1 ∑ C i=1 d(u i , U) ∑ C i=1 ∑ x j ∈CL i Ed x j , u i , ( 6
)
where N S is the number of samples, C the number of clusters, Ed the Euclidian distance between two elements, u i the i th centroid, U the center of gravity of the whole dataset, CL i the i th cluster, x j is the j th element of the dataset. I CH is the maximum value obtained for a given maximum number of clusters C. Time series of water levels were then obtained computing the median of the consecutive acquisitions from river and wetland classes. In the case of long continuous classes, time series were computed on a maximum distance of 5 km.

The silhouette coefficient provides an assessment of partitioning validity and can be used for determining the optimal number of classes [START_REF] Rousseeuw | Silhouettes: A graphical aid to the interpretation and validation of cluster analysis[END_REF].

The validity of the clustering is measured by computing the similarity or cohesion of each cluster member and its separation to the other clusters. The silhouette coefficient for the i th cluster member is defined as follows:

s(i) = b(i) -a(i) max{a(i), b(i)} , (7)
where a(i) is the average Ed between the i th member in cluster A and the other members of this cluster; b(i) is the average Ed between member i and the members in its second closest cluster B.

The silhouette ranges between -1 and 1. If s(i) close to 1, the member i is well-assigned to cluster A. If s(i) is close to 0, the object can be in-between of A and B. If s(i) is close to -1, the object is badly assigned to its cluster and on average closer to members of B [START_REF] Rousseeuw | Silhouettes: A graphical aid to the interpretation and validation of cluster analysis[END_REF].

Automatic Generation of Time Series of Water Levels

Based on the results of the unsupervised classification, time series of water levels are generated from RA data over the classes corresponding to open water and inundated areas. A time series is generated if there are at least five non-flagged altimetry heights within a distance of 3 km non-separated between each other by more than 1 km. The minimum distance between two different VSs is 3 km.

Validation

Two types of evaluation of the classification of altimetry data are performed. Global validation is performed against the land cover map of the Cuvette Centrale from [START_REF] Betbeder | Mapping of central africa forested wetlands using remote sensing[END_REF] considered as the reference. The confusion matrix is estimated when comparing the classes derived from the k-means clustering technique to the classes of land types.

Local validation parameters (bias, RMSE, R) were also performed comparing the time series of water level automatically produced using the results from the unsupervised to time series of water levels from Hydroweb.

Surface Water Volume Estimates

Maps of surface water levels were obtained combining inundation extent derived from the land cover map of the Cuvette Centrale [START_REF] Betbeder | Mapping of central africa forested wetlands using remote sensing[END_REF] with radar altimetry-based water levels from ENVISAT and Jason-2. The methodology defined in [START_REF] Frappart | Floodplain water storage in the Negro River basin estimated from microwave remote sensing of inundation area and water levels[END_REF][START_REF] Frappart | Surface freshwater storage and dynamics in the Amazon basin during the 2005 exceptional drought[END_REF] and applied to different types of satellite images, such as SAR [START_REF] Frappart | Floodplain water storage in the Negro River basin estimated from microwave remote sensing of inundation area and water levels[END_REF], multispectral [START_REF] Normandin | Quantification of surface water volume changes in the Mackenzie Delta using satellite multi-mission data[END_REF] and passive microwaves [START_REF] Becker | Satellite-based estimates of surface water dynamics in the Congo River Basin[END_REF], was used in this study. It consists in interpolating altimetry-based water stages over the inundated areas using an inverse distance weighting (IDW) technique. The wetlands inundation map was used to simulate the inundation extent over 2006-2010 as the PALSAR images used to achieve this inundation map do not allow to reach a monthly temporal resolution. During a low water period (February to April and July-August), water extent is limited to open water and permanent inundation classes (numbered 1 and 2 in Section 3.2). For the first high water period and the rising and decline of the second and larger high water period (May-June and September-October), seasonally flooded areas during the two flood pulses (class 3) are also considered, and for the peak of high water period (November-January), all classes related to water (i.e., the previously mentioned and the seasonally flooded areas during a short time corresponding to the maximum of the largest flood pulse, class 4) are merged. A map of minimum water levels recorded during the whole common observation period between SAR and altimetry was computed using a hypsometric approach to take into account the difference of altitude between the river and the floodplain (see [START_REF] Frappart | Surface freshwater storage and dynamics in the Amazon basin during the 2005 exceptional drought[END_REF] for more details on the method).

Results

Automatic Generation of Time Series of Water Levels

The optimal number of classes of backscattering coefficients was chosen by analyzing the Calinski-Harabasz criterion and the silhouette index from the number of classes ranging from 2 to 10 for both ENVISAT and Jason-2 data. The evolution of these indices as a function of the number of class is presented in Figure 2 for ENVISAT and Jason-2 backscattering. As silhouette coefficients are far from 1 for both ENVISAT (maximum value of 0.58) and Jason-2 (maximum of 0.51) backscattering coefficients, they cannot be used for determining the optimal number of classes from altimetry backscattering coefficients. On the contrary, I CH presents a clear peak for both ENVISAT and Jason-2 for 4 and 3 classes respectively (Figure 2). Nevertheless, very similar values of I CH (varying from less than 5% from the maximum value) can be observed for a number of classes ranging from 3 to 5 for both altimeters. In the following, all the results are be presented using the number of classes varying from 3 to 5.

Unsupervised Classification Results

K-means unsupervised classification was applied to ENVISAT and Jason-2 backscattering at Ku-band considering 3 to 5 classes. Figures 3 and4 present the seasonal distribution of backscattering at Ku-band from ENVISAT and Jason-2 respectively for class numbers ranging from 3 to 5. The different classes obtained from the k-means approach exhibit a quite similar temporal behavior with a maximum occurring in December, a minimum in August, a secondary maximum in May and a secondary minimum in March. This seasonal change is in phase with the temporal variations in the water levels and discharge of the Congo River (e.g., [START_REF] Abdalla | Altimetry for the future: Building on 25 years of progress[END_REF][START_REF] Huete | Development of vegetation and soil indices for MODIS-EOS[END_REF]) but in opposite phase with the backscattering at L-band from PALSAR-1 images [START_REF] Abdalla | Altimetry for the future: Building on 25 years of progress[END_REF]. Larger seasonal amplitudes are observed for the classes characterized by a higher backscattering. Better separability is observed for 3 and 4 than for 5 classes (i.e., less overlap of the std). In this latter case, the distances between the center of the consecutive classes can be lower than the sum of the standard deviations of these two classes. This lower separability affects more the classes with the lower backscattering coefficients and more the Jason-2 data than the ENVISAT data. From the results obtained in Section 5.1, confusion matrices were obtained using 3 to 5 classes from the radar altimetry backscattering datasets. They are presented in Table A1 (3 classes), Table A2 Radar altimetry backscattering classes do not have the same distribution as land type classes. Maximum values range between 0.6 and 0.65, e.g., Jason-2 class 1 data are distributed at 0.6 and 0.62 in the open water class when considering 4 and 5 classes, respectively, at 0.62 and 0.64 for short duration flooded forests (Table 2 andTable A4 respectively). These results can appear low compared with those obtained for images classification. However, as the footprint of the radar altimeter has a radius of several kilometers, it encompasses very heterogenous surfaces in terms of roughness, soil and vegetation types, moisture content and presence or not of surface water.

Previous studies analyzing radar altimetry backscattering coefficients over land surfaces showed that high values of backscattering are observed over rivers and floodplains and low values over vegetation [START_REF] Papa | Use of the Topex-Poseidon dual-frequency radar altimeter over land surfaces[END_REF][START_REF] Fatras | Spaceborne altimetry and scatterometry backscattering signatures at C-and Ku-bands over West Africa[END_REF][START_REF] Frappart | Radar altimetry backscattering signatures at Ka, Ku, C, and S bands over West Africa[END_REF][START_REF] Frappart | An ERS-2 altimetry reprocessing compatible with ENVISAT for long-term land and ice sheets studies[END_REF][START_REF] Fatras | Impact of surface soil moisture variations on radar altimetry echoes at Ku and Ka bands in semi-arid areas[END_REF].

When merging results from predominantly water and vegetation, we obtain: The spatial distribution of the backscattering classes is presented for ENVISAT and Jason-2 over the whole Cuvette Centrale in Figure 5 and a smaller area in Figure 6 to illustrate the good agreement between the unsupervised classification and the land types.

Time Series of Water Levels and Volumes

A total of 250 and 358 time series of water levels were automatically generated over the Cuvette Centrale of Congo using Jason-2 and ENVISAT data respectively for five classes of backscattering (Figure 7). Comparisons were performed between them and the closest location (below 5 km) of a Hydroweb VS. Bias, RMSE and R were estimated for 24 Jason-2 and 32 ENVISAT VSs (Figure 8). Very good agreement was found between the automatically generated time series and the Hydroweb ones with R generally higher than 0.95 (18 out of 24 for Jason-2 and 31 out of 32 for ENVISAT), RMSE lower than 0.25 m (15 out of 24 for Jason-2 VS, 22 out 32 for ENVISAT VSs) and bias lower than 0.25 m (22 out of 24 for Jason-2 VS, 26 out of 32 for ENVISAT VSs).

Surface Water Storage

Monthly maps of surface water levels and associated water volume variations were derived interpolating altimetry-based water levels on the wetland extent at a monthly temporal resolution. To avoid differences caused by intermission biases, these maps were obtained only using ENVISAT VSs. The VSs used were the ones automatically derived over either the rivers or the rivers and the floodplains. The associated surface water volume variations were computed in both cases. The annual amplitude for the largest flood period (November-December) is around four times higher (above 800 km 3 ) and the flood peak occurs 1 to 2 months later (December-January) only using VSs on rivers then using VSs on both the rivers and the floodplains (Figure 9). These differences of volume are caused by differences in spatial patterns of the water stage. Maps of the minimum height were obtained from ENVISAT. The minimum height is either the minimum of the water level during the observation period over rivers and permanently flooded pixels, or the bottom topography of floodplains over non-permanently inundated areas. In Figure 10 are presented maps based on the VSs automatically created over the rivers and the floodplains (a), or based on VSs manually created over the rivers and made available by Hydroweb (b). As it can be seen from their difference (c), the minimum height is quite similar on the center of the Cuvette Centrale, and increases up to ~30 m on the edges of the study area. Considering the mean annual amplitude, either automatically created over rivers and floodplains (d) or manually created over rivers only (e), differences up to ~10 m, positive in the north and negative in the south can be observed (f). The results of the classification of the RA backscattering coefficients show that this parameter can be used to efficiently discriminate between water and vegetation under both ENVISAT and Jason-2 ground-tracks. Contrary to SAR or multi-spectral images, this type of sensor seems to be unable to discriminate between open water and water under vegetations in the floodplains (see the results of the confusion matrices presented in Tables 1, 2 and A1-A4). This is most likely due to:

(i) The size of the illuminated area at Ku-band (several km of diameter in low resolution mode (LRM)), the scene present in the altimeter footprint is very heterogenous and encompasses both rivers and the surrounding floodplains; (ii) The RA sensor is acquiring data at nadir. As a consequence, the power backscattered by the water dominates the radar echo. RA backscattering are much stronger over rivers and floodplains than over any other types of land cover (e.g., [START_REF] Papa | Use of the Topex-Poseidon dual-frequency radar altimeter over land surfaces[END_REF][START_REF] Frappart | Radar altimetry backscattering signatures at Ka, Ku, C, and S bands over West Africa[END_REF][START_REF] Frappart | Backscattering signatures at Ka, Ku, C and S bands from low resolution radar altimetry over land[END_REF]). Even under forest canopy, water levels can be retrieved (e.g., [START_REF] Frappart | Floodplain water storage in the Negro River basin estimated from microwave remote sensing of inundation area and water levels[END_REF] and this study).

All the classes exhibit a similar temporal behavior even for the classes corresponding to non-inundated areas. The smaller seasonal cycle observed for these classes can be attributed to an increase in soil moisture during the wet periods, which a change in the dielectrical properties of the soil, and hence a rise in σ 0 which can be detected even under a dense vegetation cover [START_REF] Fatras | Spaceborne altimetry and scatterometry backscattering signatures at C-and Ku-bands over West Africa[END_REF][START_REF] Frappart | Radar altimetry backscattering signatures at Ka, Ku, C, and S bands over West Africa[END_REF][START_REF] Frappart | Backscattering signatures at Ka, Ku, C and S bands from low resolution radar altimetry over land[END_REF][START_REF] Bonnefond | The Benefits of the Ka-Band as Evidenced from the SARAL/AltiKa Altimetric Mission: Quality Assessment and Unique Characteristics of AltiKa Data[END_REF].

This good ability allowed to automatically generate time series of water levels on the RA ground-tracks; 358/250 VSs were generated on the ENVISAT/Jason-2 ground-tracks, respectively. The ones on the river were validated against manually created time series of water levels from Hydroweb due to the lack of in-situ data for validation. Very similar results were obtained when using the approach presented in this study versus the manual creation of VSs (RMSE ≤ 0.25 m in more than 80%/90% of the 32/24 cases and R ≥ 0.95 in more than 95%/75% of the same cases for ENVISAT and Jason-2, respectively). The time series of the anomaly of water levels from ENVISAT and Jason-2, averaged over all the stations located on the river and the permanently inundated areas, and over non-permanently inundated areas are presented on Figure 11. The temporal variations are consistent with the hydrological regimes of the area: a first maximum in April-May, a secondary larger one in November-December, as well as the annual amplitude around 3-4 m (e.g., [START_REF] Huete | Development of vegetation and soil indices for MODIS-EOS[END_REF][START_REF] Crétaux | Hydrological applications of satellite altimetry: Rivers, lakes, man-made reservoirs, inundated areas[END_REF][START_REF] Macqueen | Some methods for classification and analysis of multivariate observations[END_REF]) and lower over the floodplains [START_REF] Robert | Le Congo Physique[END_REF], as it can be seen when considering the average ± 1 std. The lower mean annual amplitude is observed on time series from ENVISAT (Figure 11a,c) than Jason-2 (Figure 11b,d). This is most likely due to the longer temporal revisit period of ENVISAT (35 days) than Jason-2 (10 days) and the higher number of VSs on the ENVISAT ground-tracks (358) than Jason-2 (250).

Impact on Anomaly of Surface Water Storage

The densification of the network of VSs with the definition of new VSs over the floodplains and wetlands offer a unique opportunity to better understand the dynamics of this hydrological reservoir. Two parameters are used to characterize the surface reservoir: the discharge and the storage. The combination of surface extent with water level variations provides access to the surface water storage changes (e.g., [START_REF] Frappart | Floodplain water storage in the Negro River basin estimated from microwave remote sensing of inundation area and water levels[END_REF][START_REF] Normandin | Quantification of surface water volume changes in the Mackenzie Delta using satellite multi-mission data[END_REF]). In most of the river basins, time series of water level from RA are only recorded at VSs over rivers. This study offers an opportunity to measure the impact of the inclusion of VSs over floodplains on the estimates of surface water storage variations. Even if the real time series of inundation extent is not considered here, realistic monthly variations in this parameter were used. Note that a recent study showed that the use of a climatology of surface water extent was sufficient to retrieve accurate variations in surface water storage [START_REF] Frappart | The spatio-temporal variability of groundwater storage in the Amazon River Basin[END_REF]. Our study showed that the differences in bottom topography of the floodplains (Figure 10a-c) and in the patterns of the annual amplitude of water level (Figure 10d-f) have large consequences on the amplitude of the variations in surface water storage in the Cuvette Centrale of Congo. The lack of VSs on rivers for the upstream parts in the northwest of the study area leads to an underestimate of the bottom topography from 10 to 30 m, and an overestimate of the annual amplitude above 5 m, causing an overestimate of the storage by a factor of four. Besides, spatial patterns of the annual amplitude are more consistent with the spatial distribution of the wetlands in the Cuvette Centrale. For instance, if we consider the region centered on 18 • E and from 0 to 2 • N, high annual amplitude is observed in the Ngiri and Eulonga basins (below 1 • N), a close to zero annual amplitude above, in the Upper Ubangi basin when using the VSs on the floodplains in the interpolation (Figure 10d), in good agreement with the map of the probability of wetlands [START_REF] Bwangoy | Wetland mapping in the Congo Basin using optical and radar remotely sensed data and derived topographical indices[END_REF]. On the contrary, high annual amplitude is found on the whole area when only using the VSs on the rivers (Figure 10e). These differences account for the change in amplitude and temporality (Figure 9) between the study cases (VSs only over the rivers or on the rivers and the floodplains). The influence of the larger changes in water levels observed over the rivers is limited when increasing the number of VSs in the IDW both on the storage and its seasonality (e.g., the signal is dominated by the largest floodplains and the time-lag of the seasonal cycle of the Congo River tributaries).

Conclusions

Radar altimetry is a unique tool for measuring water levels over inland waterbodies, and especially over the non-gauged wetlands and floodplains. Identifying the presence of water in wetlands and floodplains (where and when) makes difficult the construction of VSs over these environments. In this study, applying a k-means unsupervised classification approach to RA backscattering coefficients allows to accurately identify the flooded (good detection over 90%) and non-flooded areas over a complex equatorial area, the Cuvette Centrale of Congo, where floodplains are mostly covered with vegetation. Contrary to other types of satellite products, RA, due to the size of its footprint (several kilometers of diameter at Ku-band), seems unable to discriminate between open water and permanently flooded areas or to provide information on the type of vegetation covering the non-flooded area. Using the classification results, time series of water levels were automatically created over the rivers and the floodplains of the Cuvette Centrale of Congo (358/250 on the ENVISAT and Jason-2 ground-tracks, respectively). The resulting time series exhibit temporal variations in good agreement with the hydrological regime of the study area in terms of temporal variations and annual amplitude. Very similar results to the manually created time series of water levels were obtained: (RMSE ≤ 0.25 m in more than 80%/90% of the 32/24 cases and R ≥ 0.95 in more than 95%/75% of the same cases for ENVISAT and Jason-2, respectively). The densification of the VSs network by including the floodplains has a strong implication on the water volume estimate.

Owing to the larger number of VSs used when considering the VSs on the floodplains, the interpolation of the water level maps is more constrained, and the spatial patterns of annual amplitudes are more consistent with the information on the wetland extent. With the availability of data from recent RA missions operating in SAR mode such as SENTINEL-3A and B, better discrimination based on the backscattering can be expected, allowing to continue the monitoring with a better accuracy of the time series of water levels over the floodplains. This approach is likely to be applied to any other extensive wetlands or floodplains and will contribute to a better estimate, in combination with satellite images, to the changes in surface water storage over land. As the high precision altimetry era started in the early/mid 1990s with the launch of Topex/Poseidon (1992) and ERS-1/2 (1991 and 1995, respectively), a long-term record of almost 30 years of surface water storage variations can soon be available. All this information will have a large importance for better understanding the hydrological processes in the floodplains and wetlands and for the comparison with data acquired by the future NASA/CNES SWOT mission, to be launched in 2022, which will the first to estimate water levels in a swath.

Figure 1 .

 1 Figure 1. (a) The Congo River Basin and its location in the climate zones of Africa from the updated Köppen-Geiger classification [43], the red rectangle encompasses the Cuvette Centrale; (b) locations of the radar altimetry ground-tracks (Envisat in red, Jason-2 in blue) inside the Cuvette Centrale. Land types come from [38].

  from the L3 Global 500 m 16-Day moderate resolution imaging spectroradiometer (MODIS) onboard the NASA Terra satellite product (MOD13A1) from 2001 to 2009; (ii) From a series of six PALSAR-1 images at 100 m of spatial resolution acquired at L-band, in HH polarization and at various viewing angles (from 18 • to 43 • ), in ScanSAR mode by the phased array type synthetic aperture radar (PALSAR) sensors onboard the Advanced Land Observing Satellite (ALOS) on 7 September 2007, 25 October 2008, 25 January 2009, 27 April 2009, 13 December 2009, and 15 March 2010.

  (i) Open water; (ii) Permanently flooded forests; (iii) Seasonally flooded forests during the two Congo River flood pulses and located alongside the river; (iv) Seasonally flooded forests during a short time corresponding to the maximum of the largest flood pulse, farther from the river; (v) Non-flooded forests.

Figure 2 .

 2 Figure 2. Evolution of the Calinski-Harabasz criterion (I CH in blue) and the silhouette (in red) criterion as a function of the number of classes for ENVISAT (upper panel); Jason-2 (lower panel) backscattering coefficients. The following classes are present in the classification of the Cuvette Centrale used in this study: (i) Open water; (ii) Permanently flooded forests; (iii) Seasonally flooded forests during the two Congo River flood pulses and located alongside the river; (iv) Seasonally flooded forests during a short time corresponding to the maximum of the largest flood pulse, farther from the river; (v) Non-flooded forests.

Figure 3 .

 3 Figure 3. Seasonal variations in the mean backscattering coefficient and its associated deviation at Ku-band from ENVISAT (2003-2010) for each class when considering (a) 3, (b) 4, (c) 5 classes.

Figure 4 .

 4 Figure 4. Seasonal variations in the mean backscattering coefficient and its associated deviation at Ku-band from Jason-2 (2008-2016) for each class when considering (a) 3, (b) 4, (c) 5 classes.

-

  The results 0.89, 0.91 and 0.92 for ENVISAT 0.89, 0.93, 0.98 for Jason-2 when considering the sum of the backscattering values in class 1 over open water, permanently and long duration flooded forests; -The results 0.72, 0.89, 0.93 for ENVISAT and 0.85, 0.79, 0.90 for Jason-2 when considering the sum of the backscattering values in the last class over short duration flooded and non-flooded forests.

Figure 5 .

 5 Figure 5. Spatial distribution of the backscattering classes (Ku band) from (a) ENVISAT and (b) Jason-2 in the Cuvette Centrale of Congo.

Figure 6 .

 6 Figure 6. Same as Figure 5 for a smaller area.

Figure 7 .

 7 Figure 7. Locations of the VSs automatically created on the Jason-2 (blue triangles) and ENVISAT (light blue squares) ground-tracks and of the VSs from Hydroweb (Jason-2, green triangles; and ENVISAT, red circles) used for comparison.

Figure 8 .

 8 Figure 8. Comparison between altimetry-based of water levels from Hydroweb and automatically generated (this study) in terms of (a) bias, (b) RMSE and (c) r. VS locations are represented using either a dot (Jason-2) or a triangle (ENVISAT).

Figure 9 .

 9 Figure 9. Time series of surface water storage of the Cuvette Centrale of Congo obtained combining the surface water extent and the time series of water levels from the ENVISAT VSs on the rivers (black), on the rivers and floodplains (blue). Their difference is presented in dashed green.

Figure 10 .

 10 Figure 10. Maps of minimum of height from ENVISAT based on the VSs (a) automatically created over the rivers and the floodplains or (b) manually created over the rivers and made available by Hydroweb, and (c) their difference. Maps of mean annual amplitude using the VSs (d) over the rivers and floodplains, (e) only over the rivers, and (f) their difference.

1 .

 1 Identification of Open Water on RA Ground-Tracks

Figure 11 .

 11 Figure 11. Time series of anomaly of water levels averaged over all the stations located on the river and the permanently inundated areas from (a) ENVISAT and (b) Jason-2, and over the nonpermanently inundated areas from (c) ENVISAT and (d) Jason-2.

Table 1 .

 1 (4 classes) and Table1(5 classes) for ENVISAT and TableA3(3 classes), TableA4(4 classes) and Table2(5 classes) for Jason-2. Confusion matrix between 5 classes of RA backscattering coefficients from ENVISAT and the 5 land type classes identified on the study sites.

	ENVISAT Class	Open Water	Flood. Perm.	Flood. Seas. Long	Flood. Seas. Short	Non-Flood.
	1	0.42	0.29	0.21	0.07	0.01
	2	0.19	0.26	0.36	0.17	0.02
	3	0.05	0.09	0.43	0.34	0.09
	4	0.00	0.03	0.34	0.49	0.14
	5	0.01	0.01	0.15	0.57	0.26

Table 2 .

 2 Confusion matrix between 5 classes of RA backscattering coefficients from Jason-2 and the 5 land type classes identified on the study sites.

	Jason-2 Class	Open Water	Flood. Perm.	Flood. Seas. Long	Flood. Seas. Short	Non-Flood.
	1	0.60	0.21	0.17	0.01	0.01
	2	0.21	0.29	0.31	0.17	0.02
	3	0.01	0.10	0.47	0.33	0.09
	4	0.00	0.01	0.38	0.45	0.16
	5	0.00	0.00	0.10	0.64	0.26

Table A2 .

 A2 Confusion matrix between 4 classes of RA backscattering coefficients from ENVISAT and the 5 land type classes identified on the study sites.

	1	0.33	0.30	0.26	0.10	0.01
	2	0.06	0.10	0.41	0.34	0.09
	3	0.01	0.02	0.25	0.53	0.19
	ENVISAT Class	Open Water	Flood. Perm.	Flood. Seas. Long	Flood. Seas. Short	Non-Flood.
	1	0.40	0.29	0.22	0.08	0.01
	2	0.08	0.17	0.43	0.26	0.06
	3	0.00	0.04	0.42	0.41	0.13
	4	0.00	0.00	0.11	0.64	0.25

Table A3 .

 A3 Confusion matrix between 3 classes of RA backscattering coefficients from Jason-2 and the 5 land type classes identified on the study sites.

	Jason-2 Class	Open Water	Flood. Perm.	Flood. Seas. Long	Flood. Seas. Short	Non-Flood.
	1	0.41	0.25	0.24	0.09	0.01
	2	0.01	0.07	0.47	0.35	0.10
	3	0.00	0.00	0.15	0.62	0.23

Table A4 .

 A4 Confusion matrix between 4 classes of RA backscattering coefficients from Jason-2 and the 5 land type classes identified on the study sites.

	Jason-2 Class	Open Water	Flood. Perm.	Flood. Seas. Long	Flood. Seas. Short	Non-Flood.
	1	0.48	0.26	0.19	0.06	0.01
	2	0.15	0.21	0.39	0.21	0.04
	3	0.01	0.05	0.40	0.42	0.12
	4	0.01	0.02	0.18	0.57	0.22

Acknowledgments:

We thank Sylvie Gourlet-Fleury (CIRAD) for helpful discussion during the preparation of the manuscript. We also thank three anonymous reviewers for their meaningful comments which helped us improving the quality of our manuscript.

Data Availability Statement:

Radar altimetry data can be found here: http://ctoh.legos.obs-mip. fr/applications/land_surfaces/altimetric_data/altis/altis. Time series of water levels can be found here: http://hydroweb.theia-land.fr/. Accessed on 1 April 2019.

Author Contributions: Conceptualization, F.F., J.B. and V.G.; methodology, F.F., J.D., L.B., F.S.; processing, F.F., P.Z. and R.B.; data curation, F.B., N.B., F.S.; writing-original draft preparation, F.F.; writing-review and editing, all authors; funding acquisition, F.F., J.B., V.G. and J.D. All authors have read and agreed to the published version of the manuscript.

Funding: This research was funded by CNES TOSCA grants number CASCHMIR and SWHYM. The APC was offered by MDPI.

Institutional Review Board Statement: Not applicable.

Informed Consent Statement: Not applicable.

Conflicts of Interest:

The authors declare no conflict of interest.

Appendix A Table A1. Confusion matrix between 3 classes of RA backscattering coefficients from ENVISAT and the 5 land type classes identified on the study sites.


Counterstained S. aureus biofilm GFP Swimmers Merged channel

Structural resistance

Biofilms (vs plantonic phase) provide resistance to biocides

ANR GreenSwimmers

GreenSwimmers (A.Briandet/INRAE)

• swimmers bank ( 120 swimmer strains) + images

• image analysis (A.Trubuil, MaIAGE/INRAE).

• Swim descriptor and statistics ⇒typology Linking swimmer types with underlying biofilm "foragers" (explore a lot) vs "bulldozers" (dig big pores) vs "ants" (explore again pores made by others) vs ... ⇒Do bacteria adapt their swimming to the underlying biofilm ? Are there species-specific patterns ?

Question

Goal

• Understand the link between underlyning biofilm and bacterial swimming. • Characterize and discriminate bacterial species according to their swimming strategy. |A^gradb| Let X i the position of the swimmer i, dX i = dV i dt

dV i = (γ(α(b) -||V i ||)V i + β ∇b ||∇b|| )dt + ε √ dt where α(b) := v bio + b(v mat -v bio )
defines the (b dependant) swimming speed, γ is a relaxation time, β a directional force, and ε a brownian noise.

dV i = (γ(α(b) -||V i ||)V i + β ∇b ||∇b|| )dt + ε √ dt dV i = (γ(α(b) -||V i ||)V i + β ∇b ||∇b|| )dt + ε √ dt 1 Normalized model Non-dimensioned model dv = γ (v 0 + b(v 1 -v 0 ) -v ) v v dt + β ∇b ∇b dt + η mod dt
where :

                                   γ = γV * A * v 0 = v 0 V * v 1 = v 1 V * β = β A * η mod ∼ N (0, ) = A * 1 Normalized model Non-dimensioned model dv = γ (v 0 + b(v 1 -v 0 ) -v ) v v dt + β ∇b ∇b dt + η mod dt
Acceleration equation Inference problem

dv dt := Y A = f A (θ, z, v, b, ∇b) + η mod where :      θ = (γ, v 0 , v 1 , β) the parameter vector f A (z, θ, b, ∇b) = γ (v 0 + b(v 1 -v 0 ) -v ) v v + β ∇b ∇b 1 Inference problem dv dt := Y A = f A (θ, z, v
dv dt := Y A = f A (θ, z, v, b, ∇b) + η mod
Bayesian framework : regression problem

• Define priors on θ and .

• Define the likelihood function

Y A ∼ N (f A (θ|X ), )
Implementation : Stan (via pystan) 1

1 Y A vs f A (X |θ) σ i σ(y A i )

Coupling swimmers and diffusion

Pore model

R 1 (t) R 2 (t) V(t)
• Swimmers leave tails

  , b, ∇b) + η mod Data and preprocessing Compute X := (z, v, b, ∇b) positions (assumption : no observation error). v : computation by finite difference from z b and ∇b : computed by observations and finite difference on biofilm and Y A by finite difference . 1

  • mechanistic model of bacterial swimming with interaction with underlying biofilm • bayesian inference • Deterministic part of the random walk model accounts for 40 % of the variance. • Allow to reproduce macroscopic indicators • Allow to discriminate bacteria We note c and b respectively the biocide and the biofilm bacteria concentrations. ∂ t c = div(σ(b)∇c) -d cb K + b on Ω ∇c • η = 0 on ∂Ω \ Γ, c = c in on Γ where the isotropic diffusion tensor σ(b) := σ min + b(σ max -σ min ) depends on b, and the Michaelis-Menten parameters d and K are uniform.

Counterstained S. aureus biofilm

Characterizing bacterial swimming

Inference validation

Parameter chains pair-plot S.Labarthe Modeling bacterial swimmers in biofilms 2021-08-05

Modeling bacterial swimmers in biofilms 2021-08-05 21 / 25


Keywords: gold, trace elements, platinum, palladium, cementation, cupellation, LA-ICP-MS, Late Celtic coins, Caesar. Resumen oro, elementos traza, platino, paladio, cementación, copelación, LA-ICP-MS, monedas célticas tardías, Cesar

Archaeometallurgical experiments were carried out at the platform of Melle (France) to study how the trace elements of gold, in particular platinum and palladium behave when gold is melted, cupelled and then cemented. Our tests proved that platinum and palladium are reliable tracers of ancient gold. They can therefore shed light on the provenance of the gold coined by Caesar in Rome around 46-44 BC. Did he indeed melt down the Celtic gold gained from his conquest of Gaul in 51 BC as suggested by textual sources? To give a first insight into this tricky problem LA-ICP-MS was performed on some Late Celtic coins and some of Caesar's coins in order to compare especially their platinum and palladium contents, and in fact our first results support this hypothesis.

INTRODUCTION

Platinum and palladium are impurities of ancient gold that are usually drawn on for provenance studies (Blet-Lemarquand et al. 2014a). For instance, they can be used to distinguish between different stocks of precious metals if their contents scaled to the gold concentrations are significantly different. At the same time, they can also help to confirm that gold coins were melted down to manufacture other coins if the contents of both groups of coins are consistent. Their role as gold tracers is based on the assumption that platinum and palladium do not separate from ARE PLATINUM AND PALLADIUM RELEVANT TRACERS FOR ANCIENT GOLD COINS? ARCHAEOMETALLURGICAL AND ARCHAEOMETRIC DATA

TO STUDY AN ANTIQUE NUMISMATIC PROBLEM Maryse Blet-Lemarquand*, Sylvia Nieto-Pelletier*, Florian Téreygeol**, Arnaud Suspène* gold in the course of the ancient metallurgical treatments such as smelting, melting, cupellation or cementation, because these elements are highly unreactive and have high melting points.

In fact, there has not been any systematic study to determine how the trace elements of gold behave during metallurgical processes but studies are available that reports experiments on the melting of gold ores or metal gold spiked with some trace elements in oxidising or reducing atmospheres (Raub 1995; Hauptmann et al. 1995; see a brief report in Ehser 2011). Nor has there been any study exploring how trace elements behave during cupellation of gold based alloys. The behaviour of certain elements, how ever, may be generalised from the different works carried out on the cupellation of silver alloys (McKerrell and Stevenson 1972; Pernicka and Bachmann 1983; L'Héritier et al. 2015). Unfortunately, no platinum elements were examined in these articles. Cementation experiments typically focus on the efficiency of the process for major elements and neglect the minor and trace elements (Wunderlich et al. 2014; Geçkinli et al. 2000; review in Craddock 2000c). The present paper presents new analytical and experimental data for the role of platinum and palladium as valuable gold tracers.

In addition, we explore the question of provenance of the gold minted by Caesar after he conquered Celtic Gaul in the middle of the first century BC. 1 Written sources report that he seized large quantities of precious metals during his campaigns. The validity of these sources can be tested by comparing the composition of Caesar's gold coins with that of the Late Celtic Gaul gold coins. As Late Celtic gold coins are made of gold based ternary alloys, Romans would have had to melt down and purify the Celtic alloys in order to mint their high purity gold coinage. This is why our archaeometallurgical experiments included several steps, especially cupellation and cementation, in order to replicate the chaîne opératoire that the Romans would have probably followed. Some first results of these archaeometallurgical experiments were published in 2014 in a paper which lays the foundations for the present study (Blet-Lemarquand et al. 2014b).

1

This study lies within the framework of researches devoted to the provenance of Augustean gold (Suspène et al. 2011; Blet-Lemarquand et al. 2015) and of the earliest gold coins struck in Gaul (current researches led by S. Nieto-Pelletier).

ARCHAEOMETALLURGICAL EXPERIMENTS AND LABORATORY ANALYSES OF THE EXPERIMENTAL ARTEFACTS

The archaeometallurgical treatments were conceived to resemble the melting down of Late Celtic gold coins and the refining of their gold. 2 The first step was to manufacture a gold-silvercopper alloy having the same major element composition as the gold coins struck by the Celtic tribe of the Arverni in the middle of the first century BC (Tab. 4) and that contains substantial amounts of platinum and palladium. The ternary alloy was then cupelled in order to eliminate copper and other base metals from the alloy. The gold-silver alloy obtained after cupellation was laminated and finally cemented to separate gold from silver. 3 Cementation alone should be sufficient to part gold from silver and copper. For this, the gold based alloy needs to be beaten into thin foils to expose the maximum surface area. As it is easier to beat a gold-silver alloy than a gold based ternary alloy a two stage process consisting of cupellation and then cementation was preferred.

A 19 th -century French gold coin4 was selected for the experiments because previous LA-ICP-MS analysis showed it contained about 1,400 ppm platinum and 340 ppm palladium (Tab. 1).5 This 6.40 g coin was melted down with 3.12 g of pure silver and 2.07 g of pure copper to manufacture a ternary alloy of about 50 % gold, 30 % silver and 20 % copper (Fig. 1A).

The first purification process was cupellation that requires that lead is added in the right proportions to the ternary alloy. A block of lead coming from a medieval site was selected because previous LA-ICP-MS analysis showed that this lead would not contaminate our samples with platinum or palladium. For the gold-silver-copper alloy to be cupelled successfully, it requires about 10 times its weight of lead to be added to the melt as advised in published tables (Hervé, 1839: 26). Four buttons were cast that were then put into cupels made from bone ashes. The cupellations were carried out in a muffle furnace under a blast of air with temperatures varying between 950 °C and 1000 °C and four white gold buttons were obtained (Fig. 1B). The process was interrupted in the course of one of the operations because the temperature had dropped. Lead was added to the button and cupellation resumed in a new cupel.

The cementation experiments were inspired by Theophilus' recipes dating back to the 12th century AD (Hawthorne and Smith 2014) as well as by the archaeological finds from Sardis. The excavations of the workshop of Sardis provide the first archaeological evidence of the cementation of gold no later than the middle of the 6th century BC (Ramage and Craddock 2000). The four white gold buttons obtained after cupellation were accordingly laminated to a 50 micrometre thickness. 6 About 10 cycles of metal rolling and annealing were necessary to achieve this thickness. Common salt was used as active agent in the cement. 7 Twice as much ground brick was added to play the role of inert carrier even though the use of brick is not certain in Antiquity. 8 The gold foils were put alternatively with the cement in the parting vessel -a cooking pot-and the mixture was dampened with vinegar (Fig. 1C). The pot was covered with a lid, which was in fact a scorifier that had been turn over and drilled to add a temperature sensor. The process took place in an open furnace fed by charcoal that covered the cementation pot entirely and the temperature was regulated 6

The gold foils excavated in Sardis which were subjected to cementation were obviously hammered and not laminated but for the sake of convenience we preferred to laminate our samples.

7 A review of the different recipes for cements was established by Halleux (1983). 8 Craddock 2000d: 204. through an adjustable electric blower. In a second cementation experiment the pot was sealed with a lute 9 as recommended by Theophilus in his treatise (Hawthorne and Smith 2014: 109) in order to prevent the acid vapours from escaping. In fact, this lute broke in the course of the cementation. In both experiments the temperatures were maintained below the melting point of the gold-silver alloy foils of approximately 1035 °C. The first cementation lasted eight hours with temperatures maintained between 650 °C and 700 °C and the second one was lengthened to 12 hours and the temperatures were raised to 700 °C-800 °C.

The cementation experiments kept the foils intact and enhanced their golden colour. The foils obtained from the first cementation had sometimes remains of cement adhering to their surface (Fig. 1D) whereas the ones coming from the second experiment were much cleaner.

The gold based alloys were sampled at every step of the purification process in order to study their microstructure by SEM and their composition using LA-ICP-MS and SEM-EDX. The aims were to check that the metallurgical treatments were efficient and to investigate the trace element composition. The LA-ICP-MS analyses were performed using the time resolved mode that enables to reconstruct depth profiles from the surface to the interior of the object, the laser ablation being made at one spot (see for instance Fig. 4B showing a surface depletion in silver thus an enrichment in gold in a cemented foil). The calculation of the concentrations only considers the signals detected when the interior of the experimental sample is reached by the laser. The obtained results are summarized in Tab. 1. 9

The lute was made of a mixture of clay and horse dung. The cupellation proved to be efficient because copper was to a great extend eliminated from the alloy: its content was reduced from 19.4 % to about 1.6 % or even lowered to 0.1 % (Tab. 1) in the sample which was subjected to a two stage cupellation (see above). The contents of the trace elements were divided by the gold concentration for each sample in order to determine how they behaved relative to gold (Fig. 2). The ratios of most of the trace elements were reduced by a factor of 10 to 100 (Zn, As, Sb, Sn, Pb, Fe) when the gold was melted down (comparison between gold coin and Au-Ag-Cu alloy on Fig. 2A) and these parameters were still lowered during the cupellation (comparison between Au-Ag-Cu alloy and cemented alloy on Fig. 2). In contrast, the Pt/Au and Pd/ Au ratios remained unchanged between the gold coin and the cupelled alloy (Fig. 2B). It can be noticed that the cupellation left some lead (0.2 to 0.5 %) in the alloy (Tab. 1). It can be concluded that most of the trace elements were separated from the gold during melting and cupellation whereas platinum and palladium remained with adhered to gold in the course of both operations.

The cementation led to the partial elimination of silver as can be seen in several ways. The foils look more golden once they were cemented. At high magnification using a scanning electron microscope their surfaces reveal characteristics of cementation: porosities which result from the parting process (they extend along the grain boundaries and the triple points), and deposits of silver chloride formed from the attack of the silver (Fig. 3). 10 The surface of the cemented samples is also depleted in silver compared to the foils that have not been cemented: a minimum of 8-10 % silver was obtained for the foils from the first cementation and the content of sil-10 See Craddock 2000c: 181 for explanations on the chemistry of the salt cementation process.

Tab. 1. Elemental composition of the coin and of the alloys obtained by archaeometallurgical experiments (LA-ICP-MS analysis)

. Each content is the average of the concentrations calculated for 3 to 6 micro sampling. Contents in weight % (Au, Ag, Cu) or weight ppm (other elements). AVG: average, SD: standard deviation.

ver is around 5 % for the sample from the second experiment (SEM-EDX analysis). However the cross-section of a foil led us to think that the first cementation did not consistently reach the core of the foil (Fig. 4A). In fact, the depth profile analysis of the cemented foils clearly established that the surface is enriched in gold while the composition of the core remained unchanged (Fig. 4B and Tab. 1). The second cementation was more efficient than the first in that a depletion of silver was achieved also inside the foil (decrease of 5 % see Tab. 1). It seems that the cementation had an effect on the composition of some minor elements in the interior of foils. Lead contents were lowered by a factor of ten and at the same time the binary alloys were enriched in some trace elements (Sb, As, Zn, Sn, Bi), probably as a result of contaminations from the cement (Tab. 1).

The Pt/Au and Pd/Au ratios remained the same after the two cementation experiments (Tab. 1). However, we cannot strictly deduce from these analyses carried out inside the foils that cementation has no effect on these parameters because the enrichment of gold inside the foils was at best slight. But the depth profile analyses undemiably showed that the Pt/Au and Pd/Au ratios are constant from the cemented surface to the interior of the foil (Fig. 4b). Thus, cementation has no effect to these ratios. Our archaeometallurgical experiments firmly established that platinum and palladium can be relevant tracers to track a gold bullion. Consequently, platinum and palladium contents can confirm whether or not ancient gold coins were melted down and used as bullions to strike new specimen. It should be noted that most of the other element are of less reliable when trying to trace the life cycle of gold coins.

PROVENANCE OF THE GOLD COINS MINTED BY CAESAR BETWEEN 46 AND 44 BC IN ROME

During the Roman Republic, gold had only been coined in times of emergency (like the Hannibalic War or the First Civil War, Burnett 2004: 49). Caesar is actually the first Roman minting authority who regularly produced gold coins from 46-44 BC at his own mint located somewhere in Rome (Woytek 2003: 263-271). It is thought that Caesar melted down the huge quantities of gold that he had seized during his campaigns in Gaul (see Suetonius) 11 and the precious metal he took from the state's treasury in 49 BC. Pliny the Elder reported that Caesar looted in the aerarium not only 30 million sestertii or coins but also 15,000 ingots of gold and 30,000 ingots of silver (see Pliny). 12 Several scholars have already discussed the impact of the arrival of great quantities of precious metals on the Roman economy (Castelin 1977, Nash 1987: 34 or van Heesch 2005). This paper contributes new information to this discussion in the form of analytical data of Ceasar's gold coins in comparison with Late Celtic coins. In so doing, we do not claim to solve this mystery once and for all, but rather offer an important piece of new evidence about the economic ramifications of Caesar's conquests. For this purpose, 15 gold coins (aurei) minted by Caesar and believed to be representative of Caesarian coins (particularly the Hirtius issue) were selected for analyses (Tab. 2 and Fig. 5).

The composition of the gold circulating in Gaul at the time of Caesar's campaign can be gleaned on the one hand from the composition of Late Celtic coins and on the other hand from Celtic gold objects. Unfortunately, the trace element compositions of the latter are still unknown. The advantage of coins over object is that they were manufactured in large quantities and that consequently the composition of one single specimen is expected to represent the mean composition of a set of them belonging to the same series.

The Celtic coins selected for our study possibly date back to the 1 st century BC (Aubin et al. 2011; Barrandon et al. 1994; Nieto-Pelletier et al. 2011) and are attributed to different tribes: Cenomani, Arverni and Osismii (Tab. 3 and Fig. 5). The concentrations of more than 50 coins analysed using LA-ICP-MS were gathered for our study. The protocol of analysis of the coins resembles that applied to the gold foils and uses the depth profile mode to determine the composition beyond the surface that may be enriched in gold as a result of silver and copper depletion 13 and that may be «contaminated» in other elements. The laser diameter was lowered to 80 micrometres so as the craters left on the coin can not be seen with the naked eyes. Three different micro-samplings were carried out on each coin and the calculated contents are the results of the means. The three series of Late Celtic coins are made of different gold-silver-copper alloys with gold contents ranging from about 50 % to 20 % depending on the issuing authorities (Tab. 4). In 13 Three patterns of depth heterogeneities were obtained for the Celtic coins from the hoard of Laniscat, see Blet-Lemarquand et al. forthcoming. contrast, the Caesar's aurei were struck from purified gold that contains about 99.4 % gold (Tab. 4). If Caesar has indeed reused the highly debased Late Celtic coins to mint his coinage, purification processes must have necessarily been carried out to eliminate copper and silver probably by cupellation and cementation. One could wonder if it was worthwhile to purify such highly debased gold. The parting vessels discovered in Roman archaeological contexts in Britain give weight to the hypothesis of the recycling of the Iron Age precious metal by the Romans (Bayley 2009: 430). 14 As we have demonstrated, it is relevant to compare Pt/Au and Pd/Au ratios to discuss the provenance of gold (see above). The contents of platinum and palladium are normalized to gold and plotted on a binary graph for the Late Celtic and the Caesar's coins (Fig. 6). Most of the Late Celtic coins fall into an ellipse characterised by Pt/Au ratios ranging from 20 and 100, Pd/Au ratios between 8 and 40 and Pt/Pd ratios around 2.5. The two Arvernian gold coins have higher Pt and Pd contents in relation to Au but their Pt/Pd ratios are close to those of the other Late Celtic coins. It would be necessary to increase the analytical data on the Late Celtic coins especially the Arverni coins to confirm these trends. 14 For instance, J. Bayley identified parting crucibles among archaeological remains found in Chichester (Britain), a Roman settlement from immediately after the Claudian conquest of Britain in AD 43 (Bayley 2009). This site which hosted the 2nd Legion for up to two years is located near Selsey, an oppidum of the Atrebatic kingdom, where large quantities of coins and fragments of gold have been found (Bayley 2009: 430). Previous analysis established that these late Iron Age gold coins were made from ternary alloys containing around 40-50 % gold. Bayley concludes that «the late Iron Age precious metal was requisitioned in the period following the Claudian conquest and [...] (purified) so it could be subsumed into the precious metal pool the Romans had in circulation.» We are indebted to J. Bayley for directing our attention to these researches. 15

The 15 Roman coins show various Pt and Pd fingerprints stretching from low values which match the Cenomani and Osismii gold compositions to higher contents which are more consistent with the trend of the two Arvernian coins. However, three Caesar's coins seem to deviate from the Late Celtic gold because they have a lower Pt/Pd ratio of 1.3 (grey dotted line): they belong to the Plancus issue (RRC 475) or have COS QVINC on their reverse (RRC 481). This progression in the composition of the gold minted by Caesar is consistent with the dating of coins and could let us to think that Caesar also 15 The actual numbers of coins of these types are most likely more numerous than the numbers quoted in Crawford 1974 because other specimens were found since. One should also keep in mind that the assessment of the numbers of dies are very arbitrary and only indicate that the Hirtius issue was really very common whereas the Plancus and especially the Caesar issues were much less widespread. used -from 45 onwards -other supplies than the gold recovered from the Celtic coins. 16 In summary, the Hirtius issue, which represents the bulk of the gold minted by Caesar in Rome, appears to be consistent with the Late Celtic gold coins when the Pt and Pd contents are compared. At first sight, the hypothesis of recycling of the Late Celtic gold coinage by Caesar to mint his Roman coins cannot be ruled out however, as expected, other sources of gold were also used. The study is still ongoing and further analyses are planned in the framework of a next research program to complete a diachronic overview of the gold coined in the Western part of 16 It is known from the literary sources that Caesar stole precious metal from the state's treasury in 49 BC (see Pliny the Elder). But we do not know where this gold accumulated in the aerarium came from. Gold previously sent by Caesar to Rome during 58 and 49? More ancient Celtic gold? Gold coming from outside? This question remains open. the Mediterranean basin between the end of the 4th century BC and the 1st century AD 17 . This data will give us arguments to hypothesize which other gold supplies could be available for Caesar.

CONCLUSION

Archaeometallurgical experiments provided firm evidence that the Pt/Au and Pd/Au ratios remain unchanged when the gold based alloy is melted down, cupelled and then cemented. We were not able to demonstrate that cementation keeps these ratios consistent, but given the moderate temperature of the cementation process we expect this to be the case as platinum and palladium do not react with the cement at such low temperatures. Hence, platinum and palladium can serve as relevant tracers to study the provenance of gold. All other trace elements of gold, in contrast, do not survive the melting and refining of gold.

This property of platinum and palladium was applied to explore the provenance of the numerous gold coins minted by Caesar in Rome around 46-44 BC. The elemental composition of Caesar's gold coins were accordingly compared to some Late Celtic coins and at first sight their platinum and palladium chemical fingerprints partially match each other. The analytical results are consistent with the written sources. These are however preliminary results which need to be confirmed and qualified by further analyses especially on Late Celtic coins and on Celtic jewels within the framework of future research programs. The uncertainties on these ratios are mainly due to variations in the Pt and Pd contents at the scale of the micro-sampling. In fact the relative standard deviation is around 30 % for the three ablations carried out on the Caesar's coins. Despite these high fluctuations the Pt/Pd ratio remains steady for each coin. This has already been noticed for other antique gold coins (Blet-Lemarquand et al. forthcoming).

ACKNOWLEDGEMENTS

Fig. 1 .

 1 Fig. 1. The different steps of the archaeometallurgical experiments (from Blet-Lemarquand et al. 2014). From left to right, A: manufacture of the gold-silver-copper alloy; B: gold-silver alloy button obtained after cupellation; C: laminated foils made of gold-silver alloy in the cementation pot before cementation takes place; D: foils after cementation.

Fig. 2 .

 2 Fig. 2. Behaviour of the trace elements during melting and cupellation.

Fig. 3 .

 3 Fig. 3. Cemented foil after the second experiment. SEM image obtained with back scattering electron detector. Surface of the Au-Ag alloy with remains of silver chloride.

Fig. 5 .

 5 Fig. 5. Pictures of some of the coins analysed for this study. 1. Osismii; Caesar's aurei: 2. RRC 466/1 (A. Hirtius); 3. RRC 475/1b (L. Plancus); 4 RRC 481/1 (COS QVINC). Sources: IRAMAT-CEB A. Arles and G. Sarah (Celtic coin) and gallica.bnf.fr (Roman coins).

Fig. 4 .

 4 Fig. 4. Cemented foils after the first experiment, left to right. 4A: Cross-section of 6.2.a observed using a binocular lens; 4B: Depth profile analysis of 6.1.b (LA-ICP-MS analysis).

Tab. 2 .

 2 Types of the Caesar's gold coins selected for LA-ICP-MS analysis. All these coins are kept at the Bibliothèque nationale de France (BnF), département des Monnaies, médailles et antiques. RRC = Crawford 1974. The dating of coins, the number of specimen known and of dies were quoted from Roman Republican Coinage (Crawford 1974: 93-4; 888). 15 Tab. 3. Late Celtic gold coins selected for LA-ICP-MS analysis. These coins are kept at the BnF or at the Musée archéologique du Mans (hoard of Les Sablons) or at the Direction Régionale des Affaires Culturelles de Bretagne, Service Régional de l'Archéologie (hoard of Laniscat). Tab. 4. Contents of gold, silver, and copper determined for the Late Celtic coins and for the Caesar's coins. Average (AVG) and standard deviation (SD). For the results of late Celtic coins, see Aubin et al. 2011 (Cenomani); Nieto-Pelletier 2013 (Arverni) and Nieto-Pelletier et al. 2011 (Osismii). Caesar's gold: unpublished results.

  All the archaeometallurgical experiments were conducted on the national experimental archaeology platform located in Melle (France) in the framework of a PCR «Paléométallurgies et experimentations. Recherches sur les chaines de production des métaux aux périodes anciennes» 2013-2015 directed by F. Téreygeol. We are grateful to F. Duyrat and to D. Hollard (Bibliothèque nationale de France, département des Monnaies, médailles et antiques), to Y. Menez (SRA Bretagne), to B. Mandy (SRA Pays de la Loire) and to M. Thauré (Musée archéologique du Mans), for allowing us to carry out analysis on their coinages. We thank Francesca Silenzi who made the second experiment of cementation.

Fig. 6 .

 6 Fig.6. Scatterplot of the ratios of palladium and platinum to gold for the Late Celtic coins compared with the Caesar's coins (LA-ICP-MS analysis). Cenomani's results:[START_REF] Aubin | Le dépôt monétaire des Sablons, Le Mans (Sarthe) : 152 statères gaulois en or allié[END_REF] Arverni's and Osismii's gold: unpublished results. Caesar's gold: unpublished results. 18

  18

Descriptions of the ancient gold refining methods can be found in[START_REF] Halleux | Méthodes d'essai et d'affinage des alliages aurifères dans l'Antiquité et au Moyen Age[END_REF] Craddock 2000b.

The amalgamation process was not taken into account because it is not attested for Antiquity(Craddock 2000a). It was used to separate gold from platinum in Colombia during the 18th century(Morrisson et al. 1999: 125). Hence this method certainly modifies the platinum and palladium fingerprints of gold.

 4 Louis XVIII, 20 francs, 1816, bare

head.5 Information about the LA-ICP-MS analysis of ancient gold coins can be found in[START_REF] Dussubieux | LA-ICP-MS analysis of platinum-group elements and other elements of interest in ancient gold[END_REF] for the depth profile mode developed for this method see[START_REF] Gratuze | Caractérisation des alliages monétaires à base d'or[END_REF], Blet-Lemarquand et al. 2009 and Blet-Lemarquand et al. forthcoming (for gold coins) and Sarah et al. 2007 (for silver coins).

Suetonius, Vie de Jules César, LIV, 2.

Pliny, Naturalis Historia, XXXIII, 17. Unfortunately, as pointed out by Crawford, «We have no idea how large theses bars of gold and silver taken from the aerarium were»(Crawford 1974: 639).