A Delaunoy

E Prados

P Gargallo

J.-P Pons

P Sturm

Minimizing the Multi-view Stereo Reprojection Error for Triangular Surface Meshes

This article proposes a variational multi-view stereo vision method based on meshes for recovering 3D scenes (shape and radiance) from images. Our method is based on generative models and minimizes the reprojection error (difference between the observed images and the images synthesized from the reconstruction). Our contributions are twofold. 1) For the first time, we rigorously compute the gradient of the reprojection error for non smooth surfaces defined by discrete triangular meshes. The gradient correctly takes into account the visibility changes that occur when a surface moves; this forces the contours generated by the reconstructed surface to perfectly match with the apparent contours in the input images. 2) We propose an original modification of the Lambertian model to take into account deviations from the constant brightness assumption without explicitly modelling the reflectance properties of the scene or other photometric phenomena involved by the camera model. Our method is thus able to recover the shape and the diffuse radiance of non Lambertian scenes.

Introduction

Multi-view stereo is the problem of recovering the shape of scenes using cameras. Given a set of images of a scene taken from different camera positions, the goal is to reconstruct the shape S, and optionally the appearance, of the object. Since it is the inverse problem of image rendering, this problem can be modelled in a Bayesian framework by minimizing the difference between the images of the reconstructed model and the observed ones, i.e. the reprojection error. The correct variational interpretation of the Bayesian analysis yields to the minimization of some energy functional defined on the images [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF][START_REF] Jin | Variational multiframe stereo in the presence of specular reflections[END_REF][START_REF] Pons | Modelling dynamic scenes by registering multi-view image sequences[END_REF][START_REF] Soatto | Tales of shape and radiance in multi-view stereo[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF] and requires the visibility of the surface points to be carefully accounted for. Handling visibility is one of the major difficulties in the stereo vision problem. Previous work copes with this difficulty more or less elegantly. Most often the authors approximate the visibility in some pre-processing steps which can be completely prior to the whole algorithm or else inside the iterations of the minimization process [START_REF] Birkbeck | Variational shape and reflectance estimation under changing light and viewpoints[END_REF][START_REF] Faugeras | Variational-principles, surface evolution, pdes, level set methods, and the stereo problem[END_REF][START_REF] Hernandez | Multiview photometric stereo[END_REF][START_REF] Jin | Multi-view stereo reconstruction of dense shape and complex appearance[END_REF][START_REF] Pons | Modelling dynamic scenes by registering multi-view image sequences[END_REF][START_REF] Soatto | Tales of shape and radiance in multi-view stereo[END_REF][START_REF] Yu | Shape and view independent reflectance map from multiple views[END_REF].

Only recently, some authors [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF] manage to rigorously and fully account for visibility in the optimization process. In this recent work, Yezzi and Soatto [START_REF] Yezzi | Stereoscopic segmentation[END_REF] (for convex surfaces) and Gargallo et al. [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF] (in the general case) provide the exact gradient of the reprojection error. This computation is done in the mathematical framework of continuous (smooth) surfaces as used in [START_REF] Solem | A geometric formulation of gradient descent for variational problems with moving surfaces[END_REF]. This work has shown that the correct computation of the reprojection error involves a new term in the gradient which forces the apparent contours generated by the reconstructed surface to perfectly match with the apparent contours in the observed images. This makes the use of additional energy terms like ballooning [START_REF] Vogiatzis | Multiview stereo via volumetric graph-cuts and occlusion robust photo-consistency[END_REF], visual hull [START_REF] Hernández | Silhouette and stereo fusion for 3d object modeling[END_REF][START_REF] Sinha | Multi-view reconstruction using photo-consistency and exact silhouette constraints: a maximum-flow formulation[END_REF] or contours unnecessary and completely removes the minimal surface bias present in many other problem formulations.

Our work can be then seen as an extension of the previous work [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF]. One of the limitations of the gradient computed in the mathematical framework of the continuous surfaces is that, in practice we mainly deal with discrete representations of the surface. So whichever discrete representation one uses, it may be a discrete level set representation or a mesh, one has to discretize the computed continuous gradient (i.e. [START_REF] Birkbeck | Variational shape and reflectance estimation under changing light and viewpoints[END_REF][START_REF] Duan | Shape reconstruction from 3d and 2d data using pde-based deformable surfaces[END_REF][START_REF] Zaharescu | Transformesh: a topology-adaptive mesh-based approach to surface evolution[END_REF] for meshes). Also, this unavoidably introduces errors and gives rise to significant numerical difficulties. Computing the correct gradient flows for triangulated meshes is the work of [START_REF] Eckstein | Generalized surface flows for mesh processing[END_REF][START_REF] Slabaugh | Active polyhedron: surface evolution theory applied to deformable meshes[END_REF]. Nevertheless, previous authors did not use energy functionals depending on visibility (3D image segmentation [START_REF] Slabaugh | Active polyhedron: surface evolution theory applied to deformable meshes[END_REF] and 3D shape matching [START_REF] Eckstein | Generalized surface flows for mesh processing[END_REF]). In this paper, we give the exact gradient of the reprojection error for the case of discrete surfaces represented by triangulated meshes.

Another limitation of the algorithms implemented in [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF] is due to the fact that the authors assume that the scene is perfectly Lambertian. In real life, such surfaces are rare and, therefore, multi-view stereo algorithms have to be robust to non-Lambertian reflections. To improve the robustness, a common strategy is to modify input images in order to remove specular highlights [START_REF] Zickler | Color subspaces as photometric invariants[END_REF]. However, these methods are strongly limited by the specific lighting configuration. Some authors do not consider the image pixels that potentially have specular reflection components: they treat these points as outliers [START_REF] Birkbeck | Variational shape and reflectance estimation under changing light and viewpoints[END_REF][START_REF] Hernandez | Multiview photometric stereo[END_REF]. To compensate for the loss of information they then need a lot of images. Other authors try to improve the robustness to non-Lambertian effects by directly incorporating a specular reflectance model in the mathematical formulation of the problem [START_REF] Jin | Multi-view stereo reconstruction of dense shape and complex appearance[END_REF][START_REF] Yu | Shape and view independent reflectance map from multiple views[END_REF], nevertheless until now, these formulations lead to very complex and consuming optimization processes which return rather inaccurate results. Also, this kind of formulation is necessarily limited to specific materials (having the same reflectance as the one used in the modelling). Another widespread idea is to use robust similarity measures [START_REF] Jin | Variational multiframe stereo in the presence of specular reflections[END_REF][START_REF] Pons | Modelling dynamic scenes by registering multi-view image sequences[END_REF][START_REF] Yang | Dealing with textureless regions and specular highlightsa progressive space carving scheme using a novel photo-consistency measure[END_REF]. However, these are rather complex to implement and cannot deal with strongly erroneous information as in the case of camera occlusions. So the second contribution of this paper is a modification of the Lambertian model in order to take into account the deviations from the constant brightness assumption based on a smooth image-based component. We do not explicitly model the reflectance properties of the scene (contrary to [START_REF] Jin | Multi-view stereo reconstruction of dense shape and complex appearance[END_REF][START_REF] Yu | Shape and view independent reflectance map from multiple views[END_REF]) and other complex photometric phenomena involved by the cameras. Thus the method is fully generic, simple, easy to implement and very efficient. It recovers quite robustly the shape and the diffuse radiance of non Lambertian scenes viewed from flawed cameras.

Problem Formulation

Inferring the 3D shape and photometric information (such as radiance) of a scene from images of that scene is the inverse problem of image synthesis. Also if one is able to generate images from a candidate model, then, the comparison of these generated images and the observed images (i.e. the data) provides a natural reprojection error measure. In fact, this approach corresponds exactly to the Bayesian formulation of the problem which consists in maximizing the likelihood (i.e. here the reprojection error) and the prior probability (priors on the model) instead of maximizing the posterior probability directly. The formulation in terms of energies corresponds to the minimization of the sum of the data term (the reprojection error) and a prior term.

Modelling of the Reprojection Error

In order to be able to compare the whole observed images (data) with the images generated by the model, it is crucial to define and model the background. Also, as shown by [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF], this allows us to be sure that the estimated foreground surface does not shrink to an empty set (which is the global optimum for most cost functionals used in other work) and this allows that contours generated by the recovered object perfectly match the apparent contours of the input images. Whereas most of the previous work assumes that the background is known (e.g. simply modelling it as uniformly black, or by exploiting given silhouette images), here we also estimate the background images B i : I i → R 3 , under the single assumption that these images are smooth, similarly as [START_REF] Yezzi | Stereoscopic segmentation[END_REF]. Here i corresponds to the index of the camera and I i is its image domain.

Now, let us assume that the scene surface S is Lambertian and the illumination static. Let C : S → R 3 be the radiance function that associates colors to the points on the surface. Ideally, the color I i (u) observed at pixel u of image I i should be equal to the color C(π -1 i,S (u)) of its backprojection π -1 i,S (u) onto the surface or, in the case where u / ∈ π i (S), to the color B i (u) of the same pixel on the background images (π i denoting the projection associated with camera i). Thus, the reprojection error of the surface is

E data = 1 2 i πi(S) I i (u) -C(π -1 i,S (u)) 2 du + Ii-πi(S) I i (u) -B i (u) 2 du .
(1) Finally, in order to well pose the problem, we use as a prior on S an additional smoothing area energy E RS = S dσ. We also assume that the background images are smooth by adding the term

E RB = i Ii |∇B i (u)| 2 du.

Breaking Down the Constant Brightness Assumption

Few, if any, real world scenes are truly Lambertian. In this section, we present a modification of the Lambertian model (1) to take into account deviations from the constant brightness assumption without explicitly modelling the reflectance properties of the scene. We consider predicted images that would be produced by a perfect Lambertian object, plus an additional corrective component. This additional term contains the non-Lambertian component of the radiance. Also, we estimate both the Lambertian color C and the corrective component F i : I i → R 3 . In order to constrain the problem, we assume the corrective component to be smooth by adding the energy term

E RF = i Ii |∇F i (u)| 2 du.
The surface radiance is estimated using the corrected image R i :

I i → R 3 such that R i = I i -F i .
In other words, in equation ( 1) we replace I i by R i .

Additionally to correct non-Lambertian reflectance effects, images F i can also absorb other complex photometric phenomena neglected in the camera model, which makes our method able to recover quite robustly the shape and the diffuse radiance C of non-Lambertian scenes viewed from flawed cameras. Contrary to other methods which try to be robust to specular effects by modelling the reflectance using parametric models [START_REF] Yu | Shape and view independent reflectance map from multiple views[END_REF], the algorithm we propose here is quite generic and can deal with many materials without any modification. Also, it does not need complex similarity measures such as normalized cross-correlation or mutual information as [START_REF] Pons | Modelling dynamic scenes by registering multi-view image sequences[END_REF], is simpler (at least from the implementation point of view), and naturally allows us to deal with significant degradations of some parts of images, as for example can be the case with occlusions of some cameras.

For optimizing our total energy E total = E data +λ S E RS +λ B E RB +λ F E RF we perform gradient descents alternatively with respect to F i , B i and S. The computation of the gradients with respect to F i and B i is classical since they are image-based [START_REF] Yezzi | Stereoscopic segmentation[END_REF]. For a fixed shape S and a fixed C, we have

∇E total (F ) = - i (I i -C • π -1 i,S -F i )h + λ F ∆F i ∇E total (B) = - i [(I i -B i )(1 -h) + λ B ∆B i ] , ( 2
)
where h is the characteristic function that indicates if u is covered by the projection of the surface S (h(u) = 1) or not (h(u) = 0, u being explained by the background); λ is the smoothness parameter; C(x) is computed taking the mean color of the projection in the corrected image R i where x is visible. In contrast, the computation of the gradient of E total with respect to S (more exactly of E data since E RB and E RF do not depend on S; E RS is classical [START_REF] Desbrun | Implicit fairing of irregular meshes using diffusion and curvature flow[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF]) is difficult and requires some work. Actually, the computation of the exact gradient of this kind of functional is the goal of [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF], where the computation is done in the purely theoretical framework of the continuous (smooth) surfaces as the one described by [START_REF] Solem | A geometric formulation of gradient descent for variational problems with moving surfaces[END_REF] using Gâteaux derivatives.

Here, contrary to [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF], we are going to directly compute the exact gradient of a surface represented by a mesh (of course the surface being non smooth, the results of [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF] do not apply). Similarly as [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF], we first rewrite the energy as an integral over the visible volume. Then in section 4 we will compute the exact gradient of E data with respect to the discrete surface mesh. For simplicity of notation, we are going to compute the gradient for a single camera and so we remove the dependency on i. For several cameras, the gradient will be the sum of the gradients associated with each camera.

Rewriting the Energy as an Integral over the Visible Volume

The first step to be able to optimize our energy functional with respect to the shape S is to rewrite the energy as an integral over the surface instead of the image. This change of variables involves adapting the measure on the surface [START_REF] Soatto | Tales of shape and radiance in multi-view stereo[END_REF] and in counting only the visible points [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF][START_REF] Pons | Modelling dynamic scenes by registering multi-view image sequences[END_REF][START_REF] Yezzi | Stereoscopic segmentation[END_REF]. This can be achieved by du = -x•n(x)

x 3 z ν S (x)dσ where dσ is the classical surface area measure, ν S is the visibility function (giving 1 for an x that is visible from the camera and 0 otherwise), n(x) is the unit normal vector to the surface S at x, and x z is the depth of x. Thus, by using the separation technique proposed by [START_REF] Yezzi | Stereoscopic segmentation[END_REF], the energy functional becomes (for a single image):

E data (S) = -S g R (x) x • n(x) x 3 z ν S (x)dσ + I-π(S) g B (u)du , = -S [g R (x) -g B (π(x))] x • n(x) x 3 z ν S (x)dσ + I g B (u)du , (3)
where

g R (x) is 1/2 [R(π(x))-C(x)] 2 and g B (u) is 1/2 [R(u)-B(u))] 2 .
The right-hand term of equation (3) does not depend on S, so in the following we intentionally omitted it as it does not contribute to the gradient expression (with respect to S). Hence, denoting g(x) = g R (x) -g B (π(x)) for convenience, the energy to minimize with respect to S is

E(S) = - S g(x) x • n(x) x 3 z ν S (x)dσ. (4)
Finally, by using Gauss' divergence theorem which states that for all domains Ω ∈ R 3 and for all vector fields w : R 3 → R 3 , Ω ∇ • wdx = ∂Ω w • ndσ, the error functional (4) can be rewritten as :

E(S) = - R 3 ∇g(x) • x x 3 z ν S (x) dx = - R 3 f (x) ν S (x) dx . ( 5
)
For convenience, we have denoted

f (x) = ∇g(x) • x x 3 z
. Let us note that above results hold for any piecewise continuous (smooth) surface.

Gradient of the Reprojection Error for Meshes

In this section we compute the discrete gradient of E(S) (and then of E data ) with respects to the surface shape, in the case where the surface S is represented as a triangulated mesh X. Vertices of X are denoted by x k and S is deformed by moving vertices x k according to some vectors V k . The method for computing the gradient of E(S) with respect to x k consists in computing the directional derivative of E(S[t]) for any deformation x k [t] = x 0 k + tV k and then in rewriting it as a scalar product of V k , i.e. as V k , G . The obtained vector G is called the gradient and the energy necessarily decreases when deforming the surface according to its opposite direction -G. Indeed, for [START_REF] Eckstein | Generalized surface flows for mesh processing[END_REF].

x k [t] = x 0 k -tG, we have [E • S] ′ (0) = -G, G ≤ 0, see
To compute the gradient, we have to understand the variations of our energy when S is deforming. Since our energy E(S) is an integral over the visible volume, its variations are only due to the variations of the visible volume. Also, as illustrated by Figure 1 (reduced to the 2D case for simplicity), when the vertex x k is moving according to V k we have to separate two cases: 1) when all the triangles adjacent to x k are visible, the variation of the visible volume is just the sum of the tetrahedra formed by the adjacent triangles and the moved point x k + V k (see Fig. 1, left). The corresponding gradient computation is detailed in Sec. 4.1. 2) when x k is generating occluding contours in images, i.e. when it is a horizon point, its movement affects the visibility of other points located behind it (called terminator points). So the variation of the visible volume is the sum of the first case term, plus the volume swept out by the crepuscular rays generated by the horizon movement (see Fig. 1, right). The corresponding gradient computation is detailed in Sec. 4.2.

Term Due to the Tetrahedra of the Visible Adjacent Triangles

In the sequel, we denote S j the j th triangle of the mesh. Following the variation of the energy caused by the visible adjacent triangles of x k , we have

E(S[t]) -E(S[0]) = j ̟ j,V k V ol[j,V,t] f (x) dx + B(t), (6)
where V ol[j, V, t] is the volume of the tetrahedron formed by the vertices of the visible triangle S j and the point x k + V k . The sign ̟ j,V k specifies if matter has been added to or removed from the object volume and it is equal to

N j • V k |V k |
, where N j is the outward surface normal of triangle S j . B(t) on the right of equation ( 6) is null except if the vertex x k is a horizon point; this additional part will be detailed in the next paragraph. Now, we parameterize the volume

V ol[j, V, t] by x(u, v, w) = x(u, v) + w tV k , where x(u, v) = x k + u ---→ x k x k1 + v ---→ x k x k2
parametrizes the triangle S j ; x k1 and x k2 being the other two vertices of

S j ; (u, v) ∈ T = {(u, v) | u ∈ [0, 1] and v ∈ [0, 1 -u]} and w ∈ [0, 1 -u -v]. By a change of variables, (6) becomes j ̟ j,V k T 1-u-v 0 f (x(u, v) + wtV k )|det( ---→ x k x k1 , ---→ x k x k2 , tV k )| dwdudv + B(t). (7) Let A i be the area of the triangle S j . Since ̟ j,V k |det( ---→ x k x k1 , ---→ x k x k2 , tV k )| = -2t A j V k • N j , equation (7) becomes -2t V k •   j A j N j T 1-u-v 0 f (x(u, v) + wtV k ) dwdudv   + B(t). (8)
It follows that the limit of E(S[t])-E(S[0]) t when t tends to zero is

V k •   -2 j A j N j T ∇g(x(u, v)) x(u, v) x(u, v) 3 z (1 -u -v) dudv   + lim t→0 B(t), (9)
the part in square brackets being the considered first part of the gradient.

Term Due to the Movement of the Horizon

Let H k,j be the vector such that [x k , x k + H k,j ] is the edge of the triangle S j generating the horizon. The volume corresponding to the movement of the horizon can be parameterized by the points y(u, v) of the triangle {x k , x k + H k,j , x k + v tV k } generated by the movement of the horizon. More rigorously, it can be parameterized as the set of points

x(u, v, r) = r y(u, v) where y(u, v) = x k + u H k,j + v tV k ; r corresponds to the depth of x in the view point direction; r ∈ [1, T (u,v) ].
Let us note that y(u, v) depends on t; we emphasize this dependency by denoting y t (u, v). By a change of variable, we get

E(S[t]) -E(S[0]) = A(t) + 1/2 (u,v)∈T r∈[1,T (u,v) ] f (r y t (u, v)) (H k,j ∧ tV k ) • y t (u, v) |y t (u, v)| r 2 dudvdr ,
where A(t) corresponds to the part described in the previous paragraph. It follows that the limit of E(S[t])-E(S[0]) t when t tends to zero is

lim t→0 A(t) + 1/2 (u,v)∈T [1,Tu] f (r y(u)) (H k,j ∧ V k ) • y(u) |y(u)| r 2 dudvdr , (10)
where y(u) = x k + u H k,j and T u = T (u,0) . Denoting L(u) = 1,Tu f (r y(u))r 2 dr, the right-hand part of equation ( 10) can be rewritten as

V k • 1 2 u∈[0,1] L(u) y(u) |y(u)| ∧ H k,j (1 -u)du . ( 11
)
Finally, remembering that f (x) = -∇g(x) x

x 3 z , we can explicit

L(u) = g(T (y(u))) -g(y(u)) 1 [y(u)] 3 z , (12)
where T (x) is the terminator of x. Then the complete gradient flow optimizing E data (with respect to x k ) is the term in brackets of lim

Experiments

We have implemented our algorithm using the Delaunay topology-adaptive meshes proposed by [START_REF] Pons | Delaunay deformable models: Topology-adaptive meshes based on the restricted delaunay triangulation[END_REF]. The visibility is computed using OpenGL Z-buffer. The evolution is done using a multi-resolution scheme and starting from the visual hull. Horizons are located using the changes of signs of the dot products of facet normals and viewpoint directions.

Tests on the Discrete Gradient

As in [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF], we first reconstruct three uniformly colored balls arranged on a plane (20 640×480 images), see Fig. 2. This way we ensure the importance of the horizon term as the color gradient is null over the surface except at interfaces between objects in images. With the only use of the interior term (section 4.1), the surface shrinks. By using only the horizon term (given in section 4.2), we correctly reconstruct and separate the balls, and apparent contours correctly reproject in images. Then we tested our algorithm on synthetic Lambertian data for the Stanford dragon mesh (Fig. 3) composed of 32 640×480 images. The result shows the correct reconstruction of the dragon, even though the texture is smooth and some parts in shadow are dark (See Tab. 1).

Non-Lambertian surfaces with photometric camera problems

In order to test the corrective model of Sec. 2.2, we modify some of the input images used in the previous tests to simulate photometric camera problems. 8 of the 32 input images have been changed by adding a various smooth colored gradient. We not only recover the shape, but we also correctly estimate the radiance of the object, cf. Fig. 4 and Tab. 1. The corrective images are recovered on the foreground region. As expected from Sec. 3, the recovered background images corresponds to the one of the input images on the complementary regions. In a second experiment, we add a specular reflectance model to the dragon images and modify 3 of the 32 images such that half of these images is black. This could simulate occluding objects in front of the camera. Results obtained with the perfect Lambertian model of Sec. 2.1 are much less accurate than the ones obtained with the corrective model of Sec. 2.2. Also, the corrective model allows to correctly separate and recover the diffuse and specular components of the images (see Fig. 5 and Tab. 1).

Finally, we tested our method on the classical Dino (16 640×480 images) and Temple (16 640×480 images) datasets (Fig. 6). We refer to the Middlebury benchmark website [START_REF] Seitz | A comparison and evaluation of multi-view stereo reconstruction algorithms[END_REF] for evaluation. Even if these datasets were intentionally chosen for being Lambertian, they exhibit a certain amount of specular reflection that can be interpreted in the resulting corrective component. This component also shows a greenish upper part and a magenta

Diffuse images

Diffuse images Specular images Specular images (Fig. 3)

Corrective model Lambertian model Corrective model (Fig. 4) (Fig. 5) (Fig. 5) Acc. 95% 0.241mm 0.199mm 0.410mm 0.300mm Comp. 0.5mm 98.3% 98.7% 97.2% 97.6%

Table 1: Numerical evaluation of the proposed method for the dragon sequence that shows accuracy at 95% and completeness at 0.5mm following [START_REF] Seitz | A comparison and evaluation of multi-view stereo reconstruction algorithms[END_REF]. (compared to ground truth) lower part which might be due to some photometric camera bias as it appears in all images.

Conclusion

In this paper we have presented a novel deformable surface approach for the multiview stereo problem based on a discrete representation. Contrary to previous works, during the evolution, we correctly deal with visibility changes by expressing the exact gradient of the reprojection error functional. In particular, exactly as in the continuous case [START_REF] Gargallo | Minimizing the reprojection error in surface reconstruction from images[END_REF], this forces the contour generators of the surface to appear at their correct location in the images and removes the minimal surface bias from which most previous methods suffer. In

Figure 1 :

 1 Figure 1: Geometric representation of the change of visibility when moving the mesh. Contrary to the interior (left), movements of the horizon (right) strongly affect the movement of the visible interface between visible and occluded volumes by creating a movement of the crepuscular rays.

t→0A 3 z( 1

 31 (t) in (9) obtained previously, plus the following term if x k belongs to the horizon : horizon edges H k,j 1 2 u∈[0,1] g(T (y(u)))-g(y(u)) y(u) ∧ H k,j |y(u)| [y(u)] -u)du . (13)

Figure 2 :

 2 Figure 2: The balls sequence. From left to right : 4 of 20 input images; results with the horizon term computed in section 4.2 from different viewpoints.

Figure 3 :

 3 Figure 3: Synthesized dragon sequence. From left to right : 1 of 32 diffuse input images; ground truth model; recovered shape by our algorithm;

Figure 4 :

 4 Figure 4: Synthesized dragon sequence. From left to right : 2 of 32 original images; components added for simulating the photometric camera effects; input images; estimated corrective component (scaled by 2 for visualization); estimated background images; estimated shape and radiance.

Figure 5 :

 5 Figure 5: Dragon sequence synthesized with specular reflectance model and occlusions. From left to right : 2 of 32 input images; recovered corrective images; ground truth corrective images; recovered shape without the corrective model; recovered shape with the corrective model.

Figure 6 :

 6 Figure 6: Dino and Temple sequence. From left to right: 1 of 16 input images; estimated non-Lambertian component (scaled by 8 for visualization); estimated background images (scaled by 2 for visualization); estimated radiance; estimated mesh seen from a different viewpoint.

addition, our method is robust to smooth non-Lambertian effects (such as smooth specular reflections) and photometric camera defects (such as photometric calibration errors and partial occlusions) by introducing an original and simple image-based corrective model which has proved to be effective on the tests we have made on synthetic and real images.

Research was supported by the Flamenco project ANR-06-MDCA-007.


L Bazin

B Lemieux-Dudon

G Siani

A Govin

A Landais

D Genty

E Michel

S Nomade

Construction of a Tephra-Based Multi-Archive Coherent Chronological Framework 1 for the Last 2 Deglaciation in the Mediterranean Region

Keywords: last deglaciation, Mediterranean region, multi-archive comparison, coherent chronology, tephra, vegetation changes

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Understanding the mechanisms driving climate changes at the Earth's surface has largely beneficiated from measurements of geochemical tracers from various climate archives. Marine cores have provided a first relatively precise timing of reference for the succession of glacial and 2 interglacial periods over the last millions of years [START_REF] Imbrie | Modeling the climatic response to orbital variations[END_REF]. Later on, polar ice cores were the first to clearly evidence the abrupt climate variability characterizing the last glacial period (Dansgaard-Oeschger events, Dansgaard et al., 1993). Concomitantly, marine cores from North Atlantic revealed the occurrence of strong iceberg-rafted episodes during the last glacial period, named Heinrich events, that are related to some Dansgaard-Oeschger events [START_REF] Heinrich | Origin and consequences of cyclic ice rafting in the Northeast Atlantic Ocean during the past 130,000 years[END_REF]. Since this pioneer period in paleoclimatology, many different archives from the continent and the ocean have recorded the succession of glacial and interglacial periods as well as the millennial-scale climate variability of the last glacial periods [START_REF] Martrat | Abrupt Temperature Changes in the Western Mediterranean over the Past 250 , 000 Years[END_REF][START_REF] Jouzel | Orbital and millennial Antarctic climate variability over the past 800,000 years[END_REF][START_REF] Cheng | Ice Age Terminations[END_REF].

Obtaining records from different latitudes is essential to accurately depict the regional as well as global variability associated with past climate changes. With the improvement of climate models [START_REF] Eyring | Overview of the Coupled Model Intercomparison Project Phase 6 (CMIP6) experimental design and organization[END_REF], it is becoming clear and necessary to combine records from different archives in order to further assess the robustness of model outputs and their capacity to reproduce the duration and pace of past climate changes at different spatial scales [START_REF] Kageyama | The PMIP4 contribution to CMIP6 -Part 1: Overview and overarching analysis plan[END_REF].

Efforts have been made to combine records of marine and continental origin covering the last deglaciation [START_REF] Clark | Global climate evolution during the last deglaciation[END_REF][START_REF] Shakun | Global warming preceded by increasing carbon dioxide concentrations during the last deglaciation[END_REF][START_REF] Moreno | A compilation of Western European terrestrial records 60-8kaBP: Towards an understanding of latitudinal climatic gradients[END_REF]. This period framed between 19 and 11 ka ago corresponds to the transition between the last glacial period and the current interglacial (i.e. Holocene). It was induced by changes in insolation, associated with changes in greenhouse gases (GHG) concentration, atmospheric and ocean circulation reorganisations as well as the melting of polar ice sheets. The last deglaciation is widely recorded within various paleoclimatic archives. Even if numerous radiocarbon ages constrain the sequence of events of this key period, the absolute timing of changes is not always consistent from one site to another essentially because of the spatial and temporal variability of reservoir ages, the calibration of the 14 C dates [START_REF] Siani | Mediterranean Sea surface radiocarbon reservoir age changes since the last glacial maximum[END_REF][START_REF] Reimer | Intcal13 and Marine13 Radiocarbon Age Calibration Curves 0-50,000 Years Cal Bp[END_REF], but also the sample selection (e.g. shells, wood, bulk sediment; [START_REF] Müller | The role of climate in the spread of modern humans into Europe[END_REF]. The combination of climate data compilations and model simulations generally agree with the major role of GHGs and Atlantic Meridional Overturning Circulation (AMOC) variations during the last deglaciation. However, the relative leads and lags between different regions and/or type of paleoclimate archives remain difficult to assess and therefore largely unknown. Indeed, even if efforts were made in order to calibrate 14 C dates on the same reference curves [START_REF] Clark | Global climate evolution during the last deglaciation[END_REF][START_REF] Shakun | Global warming preceded by increasing carbon dioxide concentrations during the last deglaciation[END_REF], the harmonization of chronologies between the ocean and the continent, when it exists, remains largely based on alignment assumptions. For example, [START_REF] Moreno | A compilation of Western European terrestrial records 60-8kaBP: Towards an understanding of latitudinal climatic gradients[END_REF] have combined terrestrial records covering the period 60-8 ka from lake sediments, speleothems, an ice core and pollen records from marine cores. They choose not to revise the published age models. However, most of the marine cores included in their compilation possess chronologies derived from the alignment of their surface records (planktonic foraminifera δ 18 O or Sea Surface Temperature -SST) with a Greenland reference record (δ 18 Oice of NGRIP usually).

While such an assumption can be easily verified close to Greenland, thanks to the numerous volcanic ash layers [START_REF] Davies | Tracing volcanic events in the NGRIP ice-core and synchronising North Atlantic marine records during the last glacial period[END_REF][START_REF] Austin | Tracing time in the ocean: A brief review of chronological constraints (60-8 kyr) on North Atlantic marine event-based stratigraphies[END_REF], it can be questionable down to the Iberian margin and the Mediterranean region, especially when there is no assessment of the synchronization uncertainty.

Here we propose an assessment of relative climate changes associated with the last deglaciation using a coherent chronological framework based on tephra layers that are used as anchors to synchronize records. These volcanic derived deposits are by nature independent from climatic alignment assumptions (e.g. [START_REF] Lane | Volcanic ash reveals time-transgressive abrupt climate change during the Younger Dryas[END_REF][START_REF] Bronk Ramsey | Improved age estimates for key Late Quaternary European tephra horizons in the RESET lattice[END_REF]. In order to verify the consistency of the timing of changes associated with this transition, we primarily focus our chronological efforts on the Mediterranean region (Figure 1), which presents two main advantages. First, it offers a large variety of climatic records obtained from different archives on the continent and from the sediments of the Mediterranean Sea. Second, this region is characterized by widely dispersed tephra layers found within both continental and marine sedimentary archives. Their geochemical characteristics allow us to identify common tephra horizons in different archives [START_REF] Albert | Glass geochemistry of pyroclastic deposits from the Aeolian Islands in the last 50 ka: A proximal database for tephrochronology[END_REF]. In some cases these tephra layers are even precisely and accurately dated using the 40 Ar/ 39 Ar method [START_REF] Galli | Middle to Late Pleistocene activity of the northern Matese fault system (southern Apennines, Italy)[END_REF][START_REF] Giaccio | First integrated tephrochronological record for the last ∼190 kyr from the Fucino Quaternary lacustrine succession, central Italy[END_REF]Albert et al., 2019), or by 14 C dating of the surrounding material [START_REF] Lee | Modeling the Age of the Cape Riva (Y-2) Tephra[END_REF][START_REF] Albert | Revisiting the Y-3 tephrostratigraphic marker: a new diagnostic glass geochemistry, age estimate, and details on its climatostratigraphical context[END_REF][START_REF] Bronk Ramsey | Improved age estimates for key Late Quaternary European tephra horizons in the RESET lattice[END_REF].

The consistency between chronologies of the different Mediterranean paleoclimate archives included in our study is obtained using the Datice chronological integration tool [START_REF] Lemieux-Dudon | Consistent dating for Antarctic and Greenland ice cores[END_REF]. Datice was initially developed for coherent ice core dating, and used to produce the reference chronology of ice cores (Antarctic Ice Core Chronology -AICC -2012, [START_REF] Bazin | An optimized multi-proxies, multi-site Antarctic ice and gas orbital chronology[END_REF][START_REF] Veres | The Antarctic ice core chronology (AICC2012): an optimized multi-parameter and multi-site dating approach for the last 120 thousand years[END_REF]. For the purpose of our study, we improved the Datice chronological integration tool and developed a multi-archive version allowing us to consistently build chronologies for sediment cores and speleothems in addition to ice cores.

The article is organized as follows. First we present the Datice Multi-archive tool and its new development to build chronologies for deposition-like archives. Second, we list the different sites included in this study and explain the rationale behind the construction of a tephra-based chronological framework for the last deglaciation in the Mediterranean region. The comparison of our coherent chronology with previously published age scales highlights the strength and usefulness of the Datice tool for multi-archive compilation applications. Finally, we discuss the coherency of climatic changes as recorded by different archives within the Mediterranean region and compare them with NGRIP δ 18 Oice reference record from Greenland.

Tools

We use the Datice chronological integration tool to build a multi-archive coherent chronology (https://datice-multi-archives.ipsl.fr/). Information on where and how to install Datice are given in a dedicated paragraph at the end of the article. Datice was initially developed to produce coherent chronologies for ice cores by integrating absolute and relative age constraints between several ice cores for both the ice and gas phases [START_REF] Lemieux-Dudon | Consistent dating for Antarctic and Greenland ice cores[END_REF][START_REF] Buiron | TALDICE-1 age scale of the Talos Dome deep ice core, East Antarctica[END_REF][START_REF] Bazin | An optimized multi-proxies, multi-site Antarctic ice and gas orbital chronology[END_REF][START_REF] Veres | The Antarctic ice core chronology (AICC2012): an optimized multi-parameter and multi-site dating approach for the last 120 thousand years[END_REF]. Following this work, additional developments were necessary to constrain annually resolved archives such as the Greenland ice cores back to 60 ka (GICC05 chronology, [START_REF] Andersen | The Greenland Ice Core Chronology 2005, 15 42 ka. Part 1: constructing the time scale[END_REF][START_REF] Rasmussen | A new Greenland ice core chronology for the last glacial termination[END_REF][START_REF] Vinther | A synchronized dating of three Greenland ice cores throughout the Holocene[END_REF][START_REF] Svensson | A 60 000 year Greenland stratigraphic ice core chronology[END_REF], and led to the integration of duration constraints [START_REF] Lemieux-Dudon | Implementation of counted layers for coherent ice core chronology[END_REF]. On a mathematical point of view, Datice is an inverse modelling method based on Bayesian statistics. It makes the best compromise between an initial age-depth model, and its associated uncertainty estimation, with absolute and relative age constraints obtained at various depths. The prior chronology is deduced from the integration with depth of background parameters and corresponding uncertainties, all defined by the user (Table 1). Taking into account the chronological constraints (i.e. all background parameters and age markers) of all sites, as well as the stratigraphic links defined across sites, Datice produces a single chronology that is coherent to all sites and benefits from a proper assessment and propagation of age uncertainties. Moreover, after obtaining the final chronology, uncertainties are also propagated to the parameters used to calculate the final age for the different archives (e.g. accumulation rate, deposition rate, thinning function, Lock-in depth).

For ice core dating, it is necessary to date both the ice and gas phases. In Datice, the prior ice chronology is calculated from given scenarios of snow accumulation and thinning function, while the gas chronology is deduced from the evolution with depth of the reconstructed lock-in depth (i.e. the depth at which the air is trapped within the ice lattice, typically 50-100 m under the surface, Table 1, [START_REF] Lemieux-Dudon | Consistent dating for Antarctic and Greenland ice cores[END_REF].

For the purpose of multi-archive coherent dating, we further developed a new version of Datice in order to extent the use of this tool to date paleoclimate archives such as speleothems and (terrestrial and marine) sediment cores. For these mono-phase archives, the age calculation is simpler than for ice cores. Contrary to ice, which is formed from the compaction of snow and continues to thin during burying, speleothems are not affected by post-depositional compaction effects. For sediment cores, compaction effects are generally not taken into account due to the difficulty to properly quantify or model them. Consequently, the age-depth relationship for such climate archives (i.e. speleothems and sediment cores) is directly obtained from the integration with depth, from the top to the current depth, of the deposition rate L(z) (e.g. the growth rate for speleothems or the sedimentation rate for marine or lake sediment cores in centimetres per year; equation 1).

= 1 1

As the age is parameterised as increasing with depth, the depth reference is taken at the surface of the sediment or the top of the speleothem. Similarly to the background parameters for ice cores, the error of the prior estimation of the deposition rate can be parameterised as constant or varying with depth. So far, Datice Multi-archive can only integrate continuous archives. When dealing with hiatuses (regardless of the type of archives), each section has to be integrated separately (before and after the hiatus). Users have the possibility to associate error correlation with the prior estimation of the deposition rate. Adding such an error correlation allows Datice Multi-archive to propose a final deposition rate with smoother variations than in the case of no error correlation (Table 1).

In order to obtain precise and coherent chronologies, different sets of markers can be included in Datice Multi-archive (Table 1). Some are specific to ice cores, such as gas age markers or gas stratigraphic links. However, the markers applied to the ice phase (e.g. age markers, duration constraints, stratigraphic links) can also be applied to the deposition-like archives. The uncertainty associated with all markers must be given in years and should take into account all sources of error (e.g. the identification of markers, resolution of data, measurement uncertainties, calibration). Ages deduced from 14 C measurements should be calibrated into calendar years before their integration into Datice Multi-archive. Datice Multi-archive can also add error correlation between different sets of markers if needed. The coherency of the chronology between the sites is ensured by the integration of stratigraphic links. These markers must be defined a priori between two cores and have to be associated with an uncertainty in years. The final chronology uncertainty calculated by Datice Multi-archive is then a compromise between the uncertainties of the background chronology and all age markers. Large uncertainties associated with the background parameters (e.g. sedimentation rate, accumulation rate, thinning function, Table 1) necessarily result in a large uncertainty for the prior chronology. When such large uncertainties occur, the optimisation module of Datice Multi-archive has sufficient freedom to modify accordingly the background parameters and the chronology to fit the more robust dating constraints (i.e. the age markers and the stratigraphic links). Inversely, when we are confident in our background parameters, the final chronology will remain close to its initial age-depth model.

Datice Multi-archive is a complementary tool to the already existing applications used to date paleoclimate archives. The other Bayesian tools such as OxCal, BChron or Bacon [START_REF] Haslett | A simple monotone process with application to radiocarbon-dated depth chronologies[END_REF][START_REF] Bronk Ramsey | BAYESIAN ANALYSIS OF RADIOCARBON DATES[END_REF][START_REF] Blaauw | Flexible paleoclimate age-depth models using an autoregressive gamma process[END_REF][START_REF] Parnell | A review of statistical chronology models for highresolution, proxy-based Holocene palaeoenvironmental reconstruction[END_REF], are not suited for the specific problematic of ice core dating, but they offer the possibility to calibrate 14 C ages as well as build consistent chronologies for sediment cores and/or speleothems [START_REF] Bronk Ramsey | Integrating timescales with time-transfer functions : a practical approach for an INTIMATE database[END_REF], 2015). On the other hand, IceChrono [START_REF] Parrenin | IceChrono1: A probabilistic model to compute a common and optimal chronology for several ice cores[END_REF] builds coherent chronologies for ice cores only, similarly to the previous Datice version. Thanks to our recent developments, Datice

Multi-archive is the only chronological tool allowing to build coherent chronologies for ice cores as well as other paleoclimatic archives such as speleothems or (terrestrial and marine) sediment cores.

Datice Multi-archive offers an easy way to parameterize the relative constraints between the different archives through the implementation of the stratigraphic links with the possibility of adding error correlation between markers.

Consequently, the strength of Datice Multi-archives when building chronologies for paleoclimatic archives is its capacity to reconstruct chronologies for ice cores and deposition-like archives, with the proper transfer of uncertainties. Datice Multi-archive is therefore especially useful when addressing questions about the relative timing of changes associated with major climatic events that are recorded by proxies from multiple sites and archives, such as in the case of the last deglaciation.

Material and Methods

To test this new Datice version, we choose to integrate paleoclimate archives from the Mediterranean region that cover continuously the last deglaciation. The selected archives were taken on the continent (speleothems and lake/peat sediments), and from the Mediterranean Sea, and exhibit multiple high-resolution proxies. The selected marine and lake sediment cores present numerous common and independently dated tephra layers (Figure 2). This unique set of chronological tie points is of key importance to safely discuss the relative changes between the sediment records without any assumption of climatic synchronicity, as usually implied in case of chronological alignment. Consequently, we selected one marine core (MD90-917), three continental sediment cores (Monticchio, Ohrid, Tenaghi Philippon) and two speleothems (Sofular, La Mine) spread over the Mediterranean perimeter (Figure 1, Table 2). The two selected speleothems, while being independently dated and not integrated in our tephra-based Mediterranean coherent chronological framework, bring complementary information about the expression of the last deglaciation within the Mediterranean region. They, moreover, serve as test cases to validate the developments and use of Datice Multi-archive to date deposition-like archives (see section 4.1). We summarized in Table 2 the independent radiometric constraints and relative tie-points available for each site.

Both speleothems chronologies are based on numerous 230 Th ages [START_REF] Genty | Precise dating of Dansgaard-Oeschger climate oscillations in western Europe from stalagmite data[END_REF][START_REF] Fleitmann | Timing and climatic impact of Greenland interstadials recorded in stalagmites from northern Turkey[END_REF][START_REF] Göktürk | Climate on the southern Black Sea coast during the Holocene: Implications from the Sofular Cave record[END_REF]. For Monticchio, three tephra layers corresponding to ashes from the Somma-Vesuvius volcanic field have been historically recognised and are associated with a 1-sigma uncertainty of 1 year ("*" in Figure 2). The absolute ages of the Mercato and Y2 tephra layers were obtained after combining together numerous 14 C calibrated ages from different archives [START_REF] Zanchetta | Tephrostratigraphy, chronology and climatic events of the Mediterranean basin during the Holocene: An overview[END_REF][START_REF] Lee | Modeling the Age of the Cape Riva (Y-2) Tephra[END_REF][START_REF] Bronk Ramsey | Improved age estimates for key Late Quaternary European tephra horizons in the RESET lattice[END_REF]. This method statistically reduces the uncertainty associated with one volcanic event that is well-identified and 14 C-dated in numerous paleoclimate archives, especially when no 40 Ar/ 39 Ar date is available. The remaining absolute ages correspond to 40 Ar/ 39 Ar dates obtained on inland deposits of well-known eruptions ("°" in Figure 2, [START_REF] Pappalardo | Chemical and Sr-isotopical evolution of the Phlegraean magmatic system before the Campanian Ignimbrite and the Neapolitan Yellow Tuff eruptions[END_REF][START_REF] Galli | Middle to Late Pleistocene activity of the northern Matese fault system (southern Apennines, Italy)[END_REF][START_REF] Giaccio | First integrated tephrochronological record for the last ∼190 kyr from the Fucino Quaternary lacustrine succession, central Italy[END_REF]Albert et al., 2019) Finally, the chronological coherency of all sediment cores is assured by the identification of common tephra layers, which transfer the absolute age of an independently and better dated record to a poorly dated one (Figure 2). The identifications and correlation of the tephra layers are based on the published geochemical composition of tephra and their already proposed correlations [START_REF] Siani | Mediterranean Sea surface radiocarbon reservoir age changes since the last glacial maximum[END_REF], 2004[START_REF] Wulf | Tephrochronology of the 100ka lacustrine sediment record of Lago Grande di Monticchio (southern Italy)[END_REF], 2008, 2018[START_REF] Tomlinson | Geochemistry of the Phlegraean Fields (Italy) proximal sources for major Mediterranean tephras: Implications for the dispersal of Plinian and co-ignimbritic components of explosive eruptions[END_REF], 2014;[START_REF] Leicher | First tephrostratigraphic results of the DEEP site record from Lake Ohrid (Macedonia and Albania)[END_REF]Albert et al., 2019). We further confirmed theses assignments using the RESET database [START_REF] Bronk Ramsey | Improved age estimates for key Late Quaternary European tephra horizons in the RESET lattice[END_REF].

We assigned an uncertainty of 1-year to these volcanic stratigraphic links. Tephra deposition after an eruption can be considered instantaneous (i.e. <1 year). While we may encounter some bioturbation perturbation in the top-most centimetres of the sediment locally, the chronological error associated to this parameter remains difficult to model and quantify [START_REF] Bard | Reconstruction of the last deglaciation: deconvolved records of ?18O profiles, micropaleontological variations and accelerator mass spectrometric14C dating[END_REF][START_REF] Manighetti | Chronology for climate change: Developing age models for the biogeochemical ocean flux study cores[END_REF][START_REF] Carey | Influence of convective sedimentation on the formation of widespread tephra fall layers in the deep sea[END_REF][START_REF] Charbit | Effects of benthic transport processes on abrupt climatic changes recorded in deep-sea sediments: A time-dependent modeling approach[END_REF][START_REF] Barsanti | Sediment reworking rates in deep sediments of the Mediterranean Sea[END_REF]. However, enlarging the uncertainty of our volcanic-based stratigraphic links up to 10 years to account for these potential perturbations do not affect significantly the final chronology (not shown). Such a test gives us confidence in the tephrabased chronology using a 1-year uncertainty associated to the tephra correlation.

To summarize, our coherent chronological context of the last deglaciation in the Mediterranean region is constrained by 53 independently dated age markers and 20 relative stratigraphic links (Table 2 and Figure 2) for the sediment cores. Additional 59 and 13 absolute age markers constrain the Sofular and La Mine speleothem chronologies respectively, which remain here independent from the tephra-based common chronology.

To complete the set of a priori parameters (Table 1), we need to propose a prior estimation of the sedimentation/deposition rates and associated uncertainties for each deposition-like archive integrated into Datice Multi-archive. Ideally these estimates have to be independent from the age constraints, which is not always achievable. For Monticchio, we deduced the prior sedimentation rate from its published varved chronology covering the last deglaciation [START_REF] Allen | Rapid environmental changes in southern Europe during the last glacial period[END_REF] and assigned a 10% uncertainty, as proposed in the original paper (Table 3). For the other Mediterranean archives, no prior estimation of the sedimentation rate independent from the absolute age markers were available. Consequently we decided to estimate a constant background sedimentation rate from the linear regression of all the age markers for each site and associated it with a large constant uncertainty (>= 50%, Table 3). Such a parameterization allows Datice Multi-archive to substantially modify the background chronology and, as a result, respect the collection of age constraints. The chronology produced by Datice then better agrees with the age markers than with the prior guess. 2.

Figure 2: Common and well-dated tephra layers in the Mediterranean region identified in the selected terrestrial and marine sediment cores (Monticchio, MD90-917, Ohrid and Tenaghi Philippon -TP, [START_REF] Siani | Mediterranean Sea surface radiocarbon reservoir age changes since the last glacial maximum[END_REF], 2004[START_REF] Wulf | Tephrochronology of the 100ka lacustrine sediment record of Lago Grande di Monticchio (southern Italy)[END_REF], 2008, 2018[START_REF] Tomlinson | Geochemistry of the Phlegraean Fields (Italy) proximal sources for major Mediterranean tephras: Implications for the dispersal of Plinian and co-ignimbritic components of explosive eruptions[END_REF], 2014;[START_REF] Bronk Ramsey | Improved age estimates for key Late Quaternary European tephra horizons in the RESET lattice[END_REF][START_REF] Leicher | First tephrostratigraphic results of the DEEP site record from Lake Ohrid (Macedonia and Albania)[END_REF]Albert et al., 2019). Tephra in red are absolutely dated with the reference year of 1950 AD (* Somma-Vesuvius historical eruption, o 40 Ar/ 39 Ar dates or ∆ average 14 C ages with 1sigma errors). For more consistency, all 40 Ar/ 39 Ar ages are recalculated using the more recent ACs-2 standard and the 40 K total decay constant independent from astronomical calibration [START_REF] Renne | Joint determination of 40K decay constants and 40Ar * /40K for the Fish Canyon sanidine standard, and improved accuracy for 40Ar/39Ar geochronology[END_REF][START_REF] Niespolo | Intercalibration and age of the Alder Creek sanidine 40Ar/39Ar standard[END_REF]. Ages can therefore slightly differ from the original publications cited here. References : a - [START_REF] Wulf | Towards a detailed distal tephrostratigraphy in the Central Mediterranean: The last 20,000 yrs record of Lago Grande di Monticchio[END_REF][START_REF] Zanchetta | Tephrostratigraphy, chronology and climatic events of the Mediterranean basin during the Holocene: An overview[END_REF], c -Galli et al., 2017, d -Pappalardo et al., 1999, e -Bronk Ramsey et al., 2015, f -Albert et al., 2019, g -Giaccio et al., 2017.

Datice chronologies of deposition-like archives

We validate Datice Multi-archive developments to build chronologies for deposition-like archives, such as speleothems and sediment cores, using the Tunisian speleothem of La Mine and the Turkish speleothem So1 of Sofular cave (Figure 1, Table 2). The published chronologies of La Mine and So1 speleothems were obtained through linear interpolation between the 13 and 59 230 Th dates for La Mine and So1, with 2σ uncertainties ranging between 0.17 % -3.03 % and 0.26% -7.53%, of ages, respectively [START_REF] Genty | Precise dating of Dansgaard-Oeschger climate oscillations in western Europe from stalagmite data[END_REF][START_REF] Genty | Timing and dynamics of the last deglaciation from European and North African δ13C stalagmite profiles-comparison with Chinese and South Hemisphere stalagmites[END_REF][START_REF] Fleitmann | Timing and climatic impact of Greenland interstadials recorded in stalagmites from northern Turkey[END_REF][START_REF] Göktürk | Climate on the southern Black Sea coast during the Holocene: Implications from the Sofular Cave record[END_REF]. The speleothem of La Mine grew continuously between 5.5 and 23.0 ka, and the So1 speleothem of Sofular grew nearly continuously over the last 21.2 ka [START_REF] Genty | Precise dating of Dansgaard-Oeschger climate oscillations in western Europe from stalagmite data[END_REF][START_REF] Genty | Timing and dynamics of the last deglaciation from European and North African δ13C stalagmite profiles-comparison with Chinese and South Hemisphere stalagmites[END_REF][START_REF] Fleitmann | Timing and climatic impact of Greenland interstadials recorded in stalagmites from northern Turkey[END_REF][START_REF] Göktürk | Climate on the southern Black Sea coast during the Holocene: Implications from the Sofular Cave record[END_REF]. Following our background parameterization, the chronologies calculated by Datice Multiarchive are then mostly constrained by the age markers with a minor influence of the prior growth rate estimations.

In figures 3 and 4, we can note that the chronologies obtained with Datice are overall in very good agreement with the published age scales of La Mine and Sofular. Compared to a simple linear regression method, which leads to constant sedimentation rates in-between markers and no uncertainty propagation (e.g. red markers on panels C of Figures 3 and4), Datice Multi-archive produces smoother variations of the sedimentation rate with the calculation of the chronological uncertainty at all points (blue curves in figures 3 and 4). The resulting Datice Multi-archive chronology slightly differs from the published age scale of La Mine and Sofular over these intervals (panels B of Figures 3 and4). For La Mine (Figure 3), the largest difference between the chronologies reaches about 150 years, and occurs when Datice Multi-archive does not strictly respect the absolute age of the nearby marker. The final chronology remains, however, in agreement with the absolute dates within the uncertainty range (i.e. 170 years at 12.6 ka). The difference between both chronologies is induced by the smoother evolution of age with depth with Datice Multi-archive than for the original chronology of La Mine. This difference occurs for the first marker encountered with a significantly larger uncertainty than all three previous age constraints. In such a case, Datice allows the prior estimation of the deposition rate to have more impact on the final chronology relatively to the previous, more precise, age markers. For Sofular (Figure 4), the difference between both chronologies ranges between +150 years/-300 years for Datice Multi-archive chronology. Similarly to La Mine, the largest differences between the published and the Datice Multi-archive chronologies are observed when our chronological integration tool does not strictly respect the age markers. The first large age difference of 300 years is observed at around 0.45 m from the top of the speleothem at a time when Datice Multi-archive does not respect the marker at ~6 ka, which is associated with a larger uncertainty than the surrounding ones (226 years at 1σ, against 69 years for the two closest markers). A strict fit to this age marker at ~6 ka would have induced a strong shift in the deposition rate. Instead, Datice Multi-archive prefers to optimize and bring in agreement the chronology with the age markers presenting smaller uncertainties, preventing an abrupt change in deposition rate not primarily documented in the background estimate. Finally, the periods where the chronology uncertainties produced by Datice Multi-archive are the largest correspond to tipping periods close to a significant change in the deposition rate of calcite. In this case, Datice Multi-archive tries to reconcile the influence of two very close markers with significantly different absolute ages.

Using Datice Multi-archive to date mono-archives provides a proper assessment of the chronology uncertainty and sedimentation rate in-between the age markers. This method is more integrative than a simple linear interpolation between the age markers, as it was implemented in the Sofular and La Mine original publications. Figure 4: Same as Figure 3 for the Sofular So1 speleothem (published chronology from [START_REF] Fleitmann | Timing and climatic impact of Greenland interstadials recorded in stalagmites from northern Turkey[END_REF][START_REF] Göktürk | Climate on the southern Black Sea coast during the Holocene: Implications from the Sofular Cave record[END_REF].

Datice Multi-archive coherent chronology for the Mediterranean region

In this section we present the results of the Mediterranean coherent chronological framework integrating all the selected sediment/peat archives (Table 2). The chronological coherency between the different sedimentary sites is independent of any assumption of climatic alignment and only based on the identification of common and independently-dated tephra layers (Figure 2). We choose to keep the speleothem chronologies independent from the Mediterranean sediment cores in order to remain consistent and avoid unverified assumptions of climatic alignments.

The new chronology obtained with Datice for MD90-917 is in generally good agreement with the previously published chronology covering the last 24 ka [START_REF] Siani | Tephrostratigraphy study for the last 18,000 14C years in a deep-sea sediment sequence for the South Adriatic[END_REF][START_REF] Siani | Late glacial to Holocene planktic foraminifera bioevents and climatic record in the South Adriatic Sea[END_REF], Figure 5). The published chronology of MD90-917 is based on 21 calibrated radiocarbon ages combined with 1 tephra layer back to 24 ka and using a simple linear regression between the markers [START_REF] Siani | Late glacial to Holocene planktic foraminifera bioevents and climatic record in the South Adriatic Sea[END_REF][START_REF] Siani | Geoscientific Instrumentation Methods and Data Systems Paleohydrology reconstruction and Holocene climate variability in the South Adriatic Sea[END_REF]. Before this period, ages were extrapolated based on the last 14 C calibrated age used in this study. Compared to the published chronology, we have revised the calibration of the 14 C ages into calendar ages using the IntCal13 calibration curve [START_REF] Reimer | Intcal13 and Marine13 Radiocarbon Age Calibration Curves 0-50,000 Years Cal Bp[END_REF]. The largest difference between the two chronologies, up to 1115 years, is observed at the proximity of the Y3 tephra. A second significant difference occur at the proximity of the Biancavilla tephra (352 years difference at ~16.85±0.27 ka), common with Monticchio, and is induced by the combination of the tephra correlation and the revision of the 14 C age marker measured just below the tephra layer. The remaining differences between the two chronologies are within the uncertainty range and can be explained by the revision of the calibration of the 14 C ages and possibly the different age calculation methods (i.e. best compromise for Datice Multi-archive vs. linear interpolation for the published chronology). The close resemblance of the two chronologies confirms the reservoir ages estimates of [START_REF] Siani | Mediterranean Sea surface radiocarbon reservoir age changes since the last glacial maximum[END_REF] for the Adriatic Sea from the late glacial to the Holocene. The final uncertainty of MD90-917 chronology is strongly modulated by the occurrence of age markers. Indeed, the largest uncertainty values are observed at periods without age makers (e.g. ~480 years between 21-30 ka and 30-40 ka, 1σ). The new Datice Multi-archive chronology for MD90-917 is updated and more robust than the published chronology of [START_REF] Siani | Late glacial to Holocene planktic foraminifera bioevents and climatic record in the South Adriatic Sea[END_REF] over the last 24 ka, especially through providing a comprehensive estimation of the chronological uncertainty and the transfer of Monticchio's absolute constraints through the numerous stratigraphic links between the two cores.

For Monticchio, the tephra-based chronology appears a bit different from the varve-counted chronology covering the last deglaciation proposed by [START_REF] Allen | Rapid environmental changes in southern Europe during the last glacial period[END_REF], which presents an uncertainty of 5-10%. The timing of the vegetation changes corresponding roughly to the onset of the Bølling Allerød can be considered consistent within the uncertainty range of both chronologies (i.e. 15.30±0.77 ka for the chronology of [START_REF] Allen | Rapid environmental changes in southern Europe during the last glacial period[END_REF] and 15.54±0.14 ka for Datice Multiarchive chronology as recorded by the temperate taxa, Figure 5). However, further back in time the varve-counted chronology and the new Datice Multi-archive age scale start to significantly differ at around 18 ka (Figure 5), at the proximity of the Verdoline tephra layer and down to the bottom of the core section. Due to the correlated tephra layers, the absolute age of Monticchio is indirectly constrained by the 14 C calibrated ages of MD90-917 from the Tufi Biancastri (GM1 -14.55±0.60 ka) back to the Pomici di base tephra (21.27±0.30 ka on the new chronology), leading to older ages for this core section. The older Datice Multi-archive chronology is further confirmed and constrained by the absolutely dated Y3 tephra layer, common for all the sediment cores included here. Similarly to MD90-917, the uncertainty of the Datice Multi-archive chronology is larger when no marker constrain the chronology, e.g. up to ~500 years between 20 ka and 30 ka and up to ~680 years back to 38 ka, but remain smaller than the 5-10% estimated uncertainty of the varved chronology.

At Ohrid, the chronology of [START_REF] Leicher | First tephrostratigraphic results of the DEEP site record from Lake Ohrid (Macedonia and Albania)[END_REF] was obtained after linear interpolation between the absolute ages of numerous tephra layers back to 637 ka. We used the same tephra layers to constrain our chronology over our period of interest (i.e. the Mercato, Y3 and Y5 tephras, Figure 2) with revised ages using harmonized 40 K decay constant and the recent ACs-2 standard for all 40 Ar/ 39 Ar dated tephra (see section 3). Both chronologies are very consistent, always in agreement within uncertainties. Moreover, similarly to the other sediment cores, the final uncertainty of the chronology tends to increase away from the age markers (Figure 5).

The original chronology of Tenaghi Philippon was obtained after the calibration of 20 14 C ages from wood or peat bulk samples following the method of [START_REF] Weninger | A 14C age calibration curve for the last 60 ka: the Greenland-Hulu U/Th timescale and its impact on understanding the Middle to Upper Paleolithic transition in Western Eurasia[END_REF]. The Datice Multiarchive chronology of Tenaghi Philippon is constrained with the revised 14 C calibrated ages proposed by Wulf et al., 2018, combined with the revised age of the Y2 tephra layer (black markers on panel D of Figure 5, [START_REF] Bronk Ramsey | Improved age estimates for key Late Quaternary European tephra horizons in the RESET lattice[END_REF]. The coherency of its chronology with the other Mediterranean sites is only assured at the proximity of the E1, Y3 and Y5 tephra layers (Figures 2 and5). When comparing both chronologies, we notice two periods with major differences: ~465 years at around 16 ka and up to 1648 years between the Y3 and Y5 tephra layers. The differences for the older period most probably originate from the influence of the tephra layers and their revised absolute ages. For the younger period, it occurs between two age markers separated by ~5 ka, where the background parameter can have an influence on the final chronology. However, when considering both chronologies uncertainties, changes remain consistent.

By comparing our tephra-based coherent chronological framework with the published chronologies of each sediment site, we present the potential of Datice Multi-archive as a tool to build coherent chronologies that are common to various paleoclimate archives. The common tephra-based chronology allows us now to discuss the relative timing of changes during the last deglaciation as recorded by different proxies from multiple climate archives in the Mediterranean region. The δ 18 Oice of NGRIP in Greenland records regional air temperature variability over the last glacialinterglacial cycle with an annual/seasonal resolution (NGRIP community Members, 2004). This record is dated back to 60 ka using annual layer counting [START_REF] Andersen | The Greenland Ice Core Chronology 2005, 15 42 ka. Part 1: constructing the time scale[END_REF][START_REF] Rasmussen | A new Greenland ice core chronology for the last glacial termination[END_REF][START_REF] Vinther | A synchronized dating of three Greenland ice cores throughout the Holocene[END_REF][START_REF] Svensson | A 60 000 year Greenland stratigraphic ice core chronology[END_REF] and is often used as a reference for comparison with new data. Compared to previous studies, we can now compare the absolute timing of changes recorded in Greenland and within the Mediterranean region independently of any alignment assumptions, thanks to our precise and coherent volcanic-based chronological framework. In figure 6, we compare the absolute timing of changes as recorded by Greenland δ 18 Oice, to the Mediterranean SST reconstruction of core MD90-917 presented on the tephra-based coherent chronology, as well as La Mine and Sofular speleothems records (panels A, B, G and H of Figure 6 respectively). Consistent ages are observed for the sharp transitions recorded by the δ 18 Oice record of NGRIP and the MD90-917 SST reconstruction characterizing the Older Dryas (OD) -Bølling Allerød (BA), BA -Younger Dryas (YD) and the YD -Holocene transitions (Figure 6). The timing of changes for MD90-917 is not significantly changed when compared to its published chronology over these periods. However, the construction of the tephra-based Datice Multi-archive chronology for MD90-917 confirms the reservoir age values proposed by [START_REF] Siani | Mediterranean Sea surface radiocarbon reservoir age changes since the last glacial maximum[END_REF] over the LGM -Holocene transition (i.e. ~820±120 years over the OD, ~390±80 years for the LGM, BA, YD and Holocene). The consistent timing of changes associated with these sharp transitions has already largely been recorded by European speleothems [START_REF] Genty | Timing and dynamics of the last deglaciation from European and North African δ13C stalagmite profiles-comparison with Chinese and South Hemisphere stalagmites[END_REF] and western Mediterranean Sea records, however relying on 14 C calibrations with constant reservoir ages and/or tuning to NGRIP records [START_REF] Cacho | Variability of the western Mediterranean Sea surface temperature during the last 25,000 years and its connection with the Northern Hemisphere climatic changes[END_REF][START_REF] Martrat | Abrupt Temperature Changes in the Western Mediterranean over the Past 250 , 000 Years[END_REF][START_REF] Jiménez-Amat | Offset timing of climate oscillations during the last two glacialinterglacial transitions connected with large-scale freshwater perturbation[END_REF]. Consequently, using an independent volcanic-based age model for the Mediterranean sediment cores, we validate the general assumption of synchronous temperature changes associated with the last deglaciation over Greenland and within the Mediterranean region, independently of any climatic alignment.

Thanks to our volcanic-based coherent chronological framework, we can discuss the relative changes of vegetation between the marine core MD90-917 [START_REF] Combourieu-Nebout | A high-resolution record of the last deglaciation in the Central Mediterranean sea: Palaeovegetation and palaeohydrological evolution[END_REF][START_REF] Lee | Modeling the Age of the Cape Riva (Y-2) Tephra[END_REF], Monticchio [START_REF] Allen | Rapid environmental changes in southern Europe during the last glacial period[END_REF], Ohrid [START_REF] Sadori | Pollen-based paleoenvironmental and paleoclimatic change at Lake Ohrid (south-eastern Europe) during the past 500 ka[END_REF] and Tenaghi Philippon [START_REF] Müller | The role of climate in the spread of modern humans into Europe[END_REF] (panels C, D, E and F of Figure 6 [START_REF] Combourieu-Nebout | A high-resolution record of the last deglaciation in the Central Mediterranean sea: Palaeovegetation and palaeohydrological evolution[END_REF]. Such individual pollen behaviour is consistent with the progressive extension of forests under wetter and possibly warmer conditions during the BA than the OD, and the continuous persistence of herbs at low altitudes. Finally, the YD -Holocene transition also indicate local differences between the pollen records (Figure 6). At Monticchio, the YD-Holocene transition corresponds to a final sharp increase in the temperate pollen record, which is now concomitant with the SST and Greenland δ 18 Oice sharp changes in our tephra-based chronology. For MD90-917 and Tenaghi Philippon, the timing of this transition remains unchanged compared to the published chronologies. Nonetheless, their pollen records present significantly different shapes in the temperate taxa records compared to Monticchio. Indeed, the YD-Holocene transition seems to be characterized by a primary increase, up to Holocene percentages in MD90-917 and halfway to Holocene amounts at Tenaghi Philippon, preceding the final shift that is synchronous with the other sites' temperature changes (Figure 6). The pollen records of Ohrid are too poorly resolved to differentiate in details the different periods characterizing the last glacial to Holocene transition. We moreover refrain from more discussion about the absolute timing of changes recorded by Tenaghi Philippon records because of the lack of stratigraphic links with the other sedimentary cores over the transition.

Depending of the local major influencing factors controlling the isotopic composition of calcite of speleothems, we can further discuss the relative changes of vegetation types based on the δ 13 Ccalc..

Due to the geological setting of La Mine speleothem, its δ 13 Ccalc. record is not significantly affected by vegetation changes over the cave because today's conditions consists in a thin soil with sparse trees and bushes [START_REF] Genty | Timing and dynamics of the last deglaciation from European and North African δ13C stalagmite profiles-comparison with Chinese and South Hemisphere stalagmites[END_REF]. Instead, the δ 13 Ccalc. at this site is more affected by changes in the soil CO2, which mostly results from the soil biogenic production, then indirectly from temperature and humidity conditions, and is controlled by the atmospheric CO2. However, it is not possible to differentiate the δ 13 Ccalc. of La Mine in terms of quantified changes in temperature, precipitation or CO2. The δ 13 Ccalc. of Sofular, on the other hand, is mostly influenced by the type and density of vegetation as well as the soil microbial activity over the cave [START_REF] Fleitmann | Timing and climatic impact of Greenland interstadials recorded in stalagmites from northern Turkey[END_REF]. The low values of δ 13 Ccalc. (~-12‰) recorded at Sofular are interpreted as a predominance of C3 plants (~trees) over the cave, while values around ~-6‰ are associated with a predominance of C4 plants (~grasses).

In order to have a broader idea of vegetation changes in the Mediterranean region during the last deglaciation, we now further compare the pollen records of the integrated sediment cores against the δ 13 Ccalc. of Sofular, which is interpreted in terms of relative dominance of C3/trees and C4/grasses plants over the cave (located at 440 m a.s.l). The Northwestern Turkey record of Sofular exhibits a predominance of trees during the Bølling Allerød, in agreement with the general trends of warmer and more humid conditions associated with this period. The δ 13 Ccalc changes are moreover simultaneous with the sharp SST transitions recorded in the Adriatic Sea. Such a behaviour remains consistent with the interpretation of the differences between the other Mediterranean pollen records in terms of altitude sensitivity. Moreover, similarly to MD90-917 and Tenaghi Philippon, the YD -Holocene transition is associated with a two-step increase in δ 13 Ccalc., and the second part of the transition is synchronous with Greenland and Mediterranean SST changes. This preceding event has already been recorded by other pollen and SST records within the western Mediterranean region [START_REF] Combourieu Nebout | Rapid climatic variability in the west Mediterranean during the last 25 000 years from high resolution pollen data[END_REF][START_REF] Jiménez-Amat | Offset timing of climate oscillations during the last two glacialinterglacial transitions connected with large-scale freshwater perturbation[END_REF].

While the principal transitions characterizing the last deglaciation in the Mediterranean region seem to be well recorded and synchronous with Greenland changes, some differences appear in terms of integrated vegetation changes between the different sites. Even if changes recorded by individual vegetation species may be synchronous with temperature changes, the predominance of taxa curves from one site to another may differ. These differences most probably correspond to altitudinal sensitivities of vegetation species and potentially different pollen transport modes. This is especially the case for sites located at high altitudes, that can remain close to shelter areas for specific species over certain periods [START_REF] Fleitmann | Timing and climatic impact of Greenland interstadials recorded in stalagmites from northern Turkey[END_REF]. Inversely, pollen records from low altitude sites are potentially more basin-wide integrated, and may be over-imprinted, compared to the higher-altitude sites, by low-altitude pollen species. Moreover, as our tephra-based coherent chronological framework evidences asynchronous changes between temperature and vegetation changes, in terms of taxa dominance, within the Mediterranean region, we recommend caution when using synthetic pollen records for climatic alignments with reference curves, such as SST or Greenland records, to derive age models for sediment cores. NGRIP community Members, 2004). All records are presented with their respective 1σ chronological uncertainty envelope. Shaded areas highlight the position of prominent events as recorded by the SST of MD90-917 (same as in Figure 5). Selected periods are indicated on top of the figure (same as Figure 5).

Conclusion

In this paper we have presented the development and validation of Datice Multi-archive in order to build common chronologies for deposition-like archives such as speleothems or sediment records.

The advantages of using Datice Multi-archive for paleoclimate studies are threefold: 1-it can simultaneously build chronologies for several sites, regardless of the type of archive (ice cores, speleothems, sediment cores), 2-it can build one single coherent and precise chronology common to all sites, and 3-it gives a proper estimation and propagation of all the chronological uncertainties.

For the first time, we have validated the use of Datice Multi-archive by focusing on the last deglaciation. We combined together records from one marine core, three lake sediment cores and two independently-dated speleothems from the Mediterranean region. We have built a coherent and precise chronology based on the identification of common tephra layers between the sediment cores. Using this chronology, independent from climatic alignment assumptions, we showed that the major climatic transitions characterizing the last deglaciation were synchronous, within uncertainties, in Greenland and within the Mediterranean region. The combination of records from sediment cores and speleothems highlights local differences in terms of the timing of predominance of vegetation taxa during the last deglaciation within the Mediterranean region. However, we point to caution when combining together records from different archives and recommend a proper assessment of the parameters affecting these records prior to multi-archive data compilation and alignment of records.

Through this first multi-archive coherent chronological context, we evidence the usefulness of Datice Multi-archive for assessing the relative timing of changes between records from different climatic archives. Datice Multi-archive can be used to build coherent chronologies from regional to global scales while focusing on specific periods for different paleoclimatic applications.

Figure 1 :

 1 Figure1: Map of all the sites discussed in this study. Speleothems are in red (A -Sofular, B -La Mine), lake/peat sediments are in green (a -Monticchio, b -Ohrid, c -Tenaghi Philippon) and the MD90-917 marine core is in blue[START_REF] Lee | Modeling the Age of the Cape Riva (Y-2) Tephra[END_REF]. All identifications are consistent with Table2.

Figure 3 :

 3 Figure3: Comparison of chronologies for La Mine. A: depth-age relationship and age makers used to constrain the chronologies. The published chronology is in red and the new chronology obtained with Datice is in blue. B: Difference between the new chronology and the published one of[START_REF] Genty | Timing and dynamics of the last deglaciation from European and North African δ13C stalagmite profiles-comparison with Chinese and South Hemisphere stalagmites[END_REF] (positive values mean that the new chronology produced by Datice is older than the published ages). C: Uncertainties of chronologies (blue for Datice, red markers for the original chronology). D: sedimentation rates. The grey curves represent the prior estimation given to Datice, the blue curves are the resulting parameters as calculated by Datice after optimization of the chronology.

Figure 5 :

 5 Figure 5: Comparison of chronologies between the published and the tephra-based coherent chronological framework (black lines = published, colored lines = tephra-based chronology with 1sigma uncertainty envelops). The positions of volcanic-based stratigraphic links are indicated by the vertical lines between the different sites. A: SST reconstruction based on planktonic foraminifera assemblages for core MD90-917 (Siani et al., 2010, 2013, this study). The top black markers indicate the position of age markers constraining its chronology. B: pollen records of Monticchio (green = temperate, orange = steppic; Allen et al., 1999). The absolute constraints for Monticchio (black markers) are present on top of the record. C: pollen records of Ohrid (green = temperate, orange = steppic ; Sadori et al., 2016). Absolute constraints are represented by the black markers on top of the records. D: pollen records of Tenaghi Philippon (green = temperate, orange = steppic; Müller et al., 2011). The position of 14 C calibrated age and absolutely dated markers of Tenaghi Philippon are indicated at the bottom of the records. Selected periods are indicated on top of the figure: YD for Younger Dryas (~13-11.5 ka), BA for Bølling Allerød (~15-13 ka), OD for Older Dryas (~16-15 ka) and LGM for Last Glacial Maximum.

Figure 6 :

 6 Figure 6: Comparison of the last deglaciation records, presented on the Mediterranean coherent chronology between 5 ka and 25 ka B1950, with NGRIP reference record δ 18 Oice (GICC05 chonology;NGRIP community Members, 2004). All records are presented with their respective 1σ chronological uncertainty envelope. Shaded areas highlight the position of prominent events as recorded by the SST of MD90-917 (same as in Figure5). Selected periods are indicated on top of the figure (same as Figure5).

Table 1 :

 1 input parameters required for Datice Multi-archive

	Background parameters		Age markers
	Ice cores	Sediment cores	Speleothems	Absolute ages
	Accumulation rate	Sedimentation rate Deposition rate	Duration constraints
	Thinning function			Stratigraphic links
	Lock-in depth			(+ gas absolute ages, gas
				stratigraphic links, delta-depth
				markers, all specific for ice cores)
	+ associated uncertainties for each parameter and each dating constraint
		(possibility of adding error correlations for paired markers)

  . As their published ages were obtained from different calibrations, we decided to homogenize all the 40 Ar/ 39 Ar ages. Ages are then by 40 Ar/ 39 Ar. They all remain in good agreement with the 14 C ages previously published.

	recalculated using the more recent ACs-2 standard of Niespolo et al., (2017) (1.1891 Ma), and based
	on the 40 K total decay constant of Renne et al., (2011) that is independent from astronomical tuning
	(Figure 2). Consequently, the corresponding ages used in this work, and presented in Figure 2, may
	slightly differ from the published ages due to this homogenization and recalibration process (see

supplementary material). Including the recent work of Albert et al., 2019, all three major eruptions of the Campi Flegrei (i.e. NYT, Y3, Y5) that are recognised at the different sites are now precisely and accurately dated

Table 2 :

 2 List of Mediterranean sites and age constraints considered here for the last deglaciation.

	ID Name	Lat. (°N)	Long. (°E)	Absolute markers	Stratigraphic links (relative)	References
	1	MD90-917 41.00 17.52	25 14 C cal. ages 2 tephras	10 to Monticchio 3 to TP 2 to Ohrid	Siani et al., 2001, 2004, 2010; this study
	a Monticchio 40.93 15.58	3 historical eruptions 4 tephras	10 to MD90-917 2 to Ohrid 1 to TP	Allen et al., 1999; Wulf et al., 2008
	b	Ohrid	41.03 20.7	3 tephras	2 to Monticchio 2 to TP	Francke et al., 2016
						2 to MD90-917
	c	Tenaghi Philippon (TP)	40.97 24.22	13 14 C cal. ages 3 tephras	2 to Ohrid 3 to MD90-917 1 to Monticchio	Müller et al., 2011; Wulf et al., 2018
							Fleitmann et al.,
	A	Sofular	41.42 31.93	59 U/Th ages		2009; Göktürk et al., 2011
	B	La Mine	35.83 9.58	13 U/Th ages		Genty et al., 2006

Table 3 :

 3 Values of constant sedimentation rates and corresponding errors used as background parameters in Datice Multi-archive.

		MD90-917 Monticchio Ohrid Tenaghi Philippon	Sofular La Mine
	Sed. rate		Varved
	(cm/a)	0.0346	chronology a	0.04	0.0338	0.008	0.001
	Sigma
		0.015	10%	0.02	0.02	0.01	0.01
	(cm/a)
	a Allen et al., 1999
	4. Results

  respectively). With our new coherent chronological framework, we now have a proper assessment of the chronological uncertainties for each site, rendering possible the inter-comparison of relative vegetation changes between Monticchio, MD90-917, Tenaghi Philippon and Ohrid. First, all records agree in indicating a predominance of temperate taxa during the Holocene, and more steppic taxa during the last glacial period. However, the timing of changes in the predominance of temperate versus steppic pollen taxa at the different sites differs significantly and is now precisely established thanks to our multi-archive coherent chronology. At Monticchio, temperate tree taxa start to increase simultaneously with the SST increase recorded in the Adriatic Sea by the MD90-917 core, and become predominant at the onset of the Bølling-Allerød period (i.e. at 14.63±0.14 ka on our coherent chronological framework). Ohrid records indicate a similar behaviour slightly later (i.e. 14.04±0.87 ka), that can be considered synchronous with Monticchio within uncertainties. On the other hand, records from MD90-917 and Tenaghi Philippon show a shift in the predominance of temperate vs. steppic taxa at ~12 ka (i.e. at 11.98±0.16 ka for MD90-917 and 12.10±0.64 ka for Tenaghi Philippon). In MD90-917 records, the BA period is characterized by a dominance of steppic taxa and also higher, but still low, amounts of temperate taxa. In contrast, Monticchio and Ohrid pollen records indicate a continuous increase in temperate pollens from the initial gradual warming recorded in the Adriatic Sea until the BA-YD sharp transition. Such regional differences could originate from the difference in altitude between the sites. While MD90-917 and Tenaghi Philippon are low-altitude sites, where pollen are mostly transported by wind and basin-wide hydrological transport, Monticchio and Ohrid sites are higher in altitudes (i.e. 656 m a.s.l. for Monticchio and 693 m a.s.l for Ohrid). They may be potentially biased by high-altitude vegetation changes and less recording vegetation changes from lower altitudes. One should note that while the temperate and steppic taxa curves present different timing of relative predominance between the four sites, the records of individual pollen species of MD90-917 (e.g. Quercus, Alnus, Betula) indicate an increase of trees in altitude on the mountainous slopes around the Adriatic Sea, while herb species remain dominant (e.g. Artemisia and Ephedra, Combourieu-

Acknowledgments

We thank the two anonymous reviewers and Anders Svensson for their comments and suggestions that have helped us to improve the manuscript during the review process. This work was supported by Labex L-IPSL, which is funded by the ANR (grant no. ANR-10-LABX-0018). This work was supported by the French national program LEFE-IMAGO/INSU. This is LSCE publication number X.

Availability:

The Datice Multi-archive software is available under a CeCILL license and can be downloaded from the INRIA Forge after registration (see https://gforge.inria.fr/account/register.php). Documentation and support for the install of Datice Multi-archive can be found on the INRIA Forge, or upon request to lucie.bazin@lsce.ipsl.fr or lemieux.benedicte@gmail.com. The git version control tool is necessary to clone the INRIA Forge repository. Compiling the executable requires a Fortran compiler with libraries BLAS-Lapack and netcdf (incidentally Openmp library is necessary to test the multi-threaded and shared memory beta version). Datice Multi-archive sources include the INRIA modulopt optimization solvers library [START_REF] Gilbert | An environment for testing solvers on heterogeneous collections of problems[END_REF]. Information on how to run Datice Multi-Archive, the different input, param and output files can be found at https://datice-multiarchives.ipsl.fr/. A visualization tool is also available when installing Datice, requiring Python (version 2.7) with modules numpy, scipy, matplotlib, wxPython and netcdf4.

Supplementary Material:

-1 file with all background parameters and age markers used to constrain the chronologies for all sites -1 file with the new Datice chronologies (depth, age, uncertainty) and data for all sites Bibliography Albert, P.G., Giaccio, B., Isaia, R., Costa, A., Niespolo, E.M., Nomade, S., Pereira, A., Renne, P.R., Hinchliffe, A., 2019. Evidence for a large-magnitude eruption from Campi Flegrei caldera ( Italy ) at 29 ka 47, 1-5. https://doi.org/10.1130/G45805.1/4684381/g45805.pdf