Marion Parquer

Guillaume Caumon

Pauline Collon

Integration of channel meander abandonment age uncertainty into a stochastic channelized system reconstruction method

Keywords: fluvial, turbidite, point bars, heterogeneity, sedimentary record, remote sensing, interpretation

published or not. The documents may come

Integration of channel meander abandonment age uncertainty into a stochastic channelized system reconstruction method

Introduction

In turbiditic or fluvial deposits, understanding the channelized system evolution is important for ancient system architecture reconstruction. However, this task is difficult because the various processes involved in the channel evolution only leave partial clues of the past geometries. Indeed, erosion and deposition by meandering channels generate architectural elements such as point bars formed by gradual channel migration (e.g., [START_REF] Gay | Evolution of cutoffs across meander necks in Powder River, Montana, USA[END_REF][START_REF] Holbrook | Expansion to contraction, translation to retreat, osculation or toggle, and degree of wobble, and the many means by which a meander may move[END_REF] and abandoned meanders due to the cutoff of meander loops (e.g., [START_REF] Lonsdale | Cut-offs at an abyssal meander south of Iceland[END_REF][START_REF] Erskine | Alluvial cutoffs as indicators of former channel conditions[END_REF].

When abandoned, the paleo-channel path (abandoned loop or avulsed channel) is generally filled with deposits of the remaining suspended sediments at the time of abandonment, and by the following successive floods occurring in the meander belt [START_REF] Toonen | Sedimentary architecture of abandoned channel fills[END_REF]. The filling is mainly composed of muddy material, both in alluvial settings [START_REF] Veeken | Seismic stratigraphy, basin analysis and reservoir characterisation[END_REF] and turbidite settings [START_REF] Issautier | Impacts of fluvial reservoir heterogeneity on connectivity: Implications in estimating geological storage capacity for CO2[END_REF]. It can hence be imaged by reflection seismic data (e.g., [START_REF] Posamentier | Seismic geomorphology and stratigraphy of depositional elements in deep-water settings[END_REF]Miall, 2014).

Thus, many observations can be made on satellite, seismic, Light Detection and Ranging (LIDAR), or Ground-Penetrating Radar (GPR) images (Fig. 1). Complete channel trajectories and abandoned meanders can be interpreted, as well as erosive patterns between them. These structures provide valuable information about the past trajectories of the main channel. In some specific cases, it is, therefore, possible to have a confident depositional history integrating all architectural elements. For example, [START_REF] Durkin | Threedimensional reconstruction of meander-belt evolution, Cretaceous McMurray Formation, Alberta Foreland Basin, Canada[END_REF] proposed an exhaustive and precise study of the point bar evolution by analyzing cross-cutting, orientation and paleo-current data. In order to determine relative or absolute ages, these geometric considerations can be completed, for example, by field work [START_REF] Fisk | Geological investigation of the alluvial valley of the lower Mississippi River[END_REF], borehole data analysis [START_REF] Hubbard | Seismic geomorphology and sedimentology of a tidally influenced river deposit, Lower Cretaceous Athabasca oil sands, Alberta, Canada[END_REF], radiocarbon analysis [START_REF] Holbrook | Stratigraphic evidence for millennial-scale temporal clustering of earthquakes on a continental-interior fault: Holocene Mississippi River floodplain deposits, New Madrid seismic zone[END_REF] or optically stimulated luminescence [START_REF] Rowland | Tie channel sedimentation rates, oxbow formation age and channel migration rate from optically stimulated luminescence (OSL) analysis of floodplain deposits[END_REF]. However, such detailed characterizations have a significant cost in terms of data acquisition and analysis. Their application to a particular system is, therefore, not always feasible, especially in subsurface settings, where direct sampling is seldom an option and where the resolution of the geophysical images, most of the time, does not permit the clear interpretation of scroll bars. Additionally, the exact time of meander abandonment or of avulsion triggering is often difficult to estimate from the sedimentary record [START_REF] Mohrig | Interpreting avulsion process from ancient alluvial sequences: Guadalope -Matarranya system (northern Spain) and Wasatch Formation (western Colorado)[END_REF][START_REF] Hajek | Significance of channel-belt clustering in alluvial basins[END_REF]. Indeed, only a relative age of an abandoned loop can be determined when cross-cutting is observed between the loops and the surrounding structures, or based on depth considerations. In contrast, isolated abandoned loops where no clear erosional pattern is observed cannot be relatively dated. This is typically a problem for some interpretation of stratal seismic slices, where depth cannot be used to infer the complete chronology of interpreted architectural elements (Figs. 1b and1c). Well data inform locally on the presence or absence of a paleo-channel facies, and occasionally about the presence of point bars and orientation of the paleo-migration [START_REF] Brekke | The use of microresistivity image logs for facies interpretations: An example in point-bar deposits of the McMurray Formation, Alberta, Canada[END_REF]. In summary, chronologic uncertainty stems from the irreversibility of abandonment, from incompleteness of the sedimentary record and from the limited quality and quantity of available data. A consistent reconstruction of the channelized system must therefore take into account all these field observations but also integrate the uncertainties that are attached to them. This calls for a stochastic simulation tool to generate 3D channel system architectures that correspond to various possible global chronologies of meander abandonments inferred from seismic or satellite images.

To address these challenges, we propose a method that uses the relative ages between architectural elements identified by cross-cuttings, and the geometry of all abandoned meander loops, which provide indication of paleo-channel trajectories. The main idea is to exploit and improve the reverse-time channel migration method (ChaRMigS (Channel Reverse-time Migration Simulation), [START_REF] Parquer | Reconstruction of channelized systems through a conditioned reverse migration method[END_REF]), both to discover the ages of the known abandoned meanders of a particular channel belt and to reconstruct possible internal geometries of the associated sedimentary deposits. Indeed, the ChaRMigS reconstruction method permits the production of equiprobable stochastic 3D geometries of a channelized system by reverse migrating the main channel path for a given number of time steps. At each time step, an abandoned meander loop can be integrated in the main channel path and then be reverse migrated. The technique presented in this paper uses a stochastic estimation of abandonment at each reverse migration time step and thus proposes possible 3D global geometries of the subsurface. Indeed, these various internal architectures are due to the various lateral and verti- cal migration simulated in the different stochastic realizations produced by ChaRMigS. Because of the inherent uncertainties about this chronology and the associated geometry of the system, the method uses a stochastic process to generate several possible realizations of paleo-channels through time.

In this paper, the starting point of the method is a full channel trajectory and anterior abandoned loops interpreted from available data. We do not use scroll bars for constraining the methodology at this stage, both for simplicity and because limited seismic resolution may prevent their identification in some cases. This approach, illustrated in Fig. 2 and Algorithm 1, exploits the relative ages that can be determined thanks to cross-cutting relationships between paleo-meanders. In the absence of an erosion relation between paleomeanders, we propose to order them by a random sampling. This sampling accounts for the spatial organization with regard to the current channel path and for cutoff statistics based on analog studies. The integration of this tool in the ChaRMigS method permits us to update the current channel path at each time step, and thus to deduce the global chronology of the abandoned loops at the same time as reconstructing the complete channelized system (Fig. 3). An application of the developed technique is presented on a seismic stratal slice of the McMurray Formation in Canada (Section 4.2).

-channel path

Relative chronology observations

As discussed above, several observations such as truncation and relative depth may indicate local relative chronologies between abandoned meanders. Constraining theses local chronologies is key to reduce the global channelized system chronologic uncertainty. In this section, we propose a semi-automatic strategy to come up with the local relative age of abandoned meanders, which uses the intersections between curves interpreted from available data, as described in Fig. 4 and in Algorithm 1.

During the evolution of a meander belt, a meander bend can erode another bend that was previously abandoned. The relation of erosion gives a relative chronology: the eroded bend is older than the eroding one. This erosion leads to the eventual truncation of the abandoned loop extremities (e.g., -compute the probability of integration (using the distance and orientation to the current channel path) end for -sample from C the number N am of abandoned meanders to be integrated at this time step t -sample the N am among the set of candidates O j according to their probabilities of integration for each group from which an abandoned meander has been integrated do for each O k not selected for integration at this time step do -decrement the local chronologies L k end for end for return one probable abandonment chronology at a given time step erosion of the neck of an abandoned neck cutoff loop). To face this issue, we simply propose a manual extrapolation of the full abandoned meander during digitization. This interpretive step could be automated in the future, but is now necessary to facilitate the integration of the abandoned meanders in the main channel path during ChaRMigS reverse migration.

Pairwise automatic intersection checks are performed between observed paleo-channel paths (Fig. 5a). The eroding/eroded abandoned meanders compose groups of abandoned loops (Fig. 5b) in which local chronologies can be deduced by visual analysis of the seismic or satellite image.

The groups of paleo-geometries indicate abandoned meanders that can be chronologically ordered through the compilation of the pairwise relative chronologies. Relative ages are assigned to the abandoned meanders in chronological order. By convention, the abandoned meanders that have been abandoned most recently have the smallest relative age. A relative age equal to zero is assigned to the youngest abandoned meander of each group. The relative age of anterior abandoned meanders is incremented within each group (Fig. 6).

An abandoned meander presenting no intersection with another observed abandoned meander is an isolated one. It constitutes a group by itself, and its local chronology value is always set to zero.

In the above semi-automatic method, local chronologies are directly de- duced by interpretation of the curve intersections inside groups. Alternatively, local relative ages could be completely provided by experts as an input to the next steps of the algorithm. However, in all cases, the global chronology between groups cannot be easily determined. During reverse migration, the knowledge (or the interpretation) of the paleo-channel trajectory provides a possible way to estimate a global chronology for abandoned loops.

Simulation of the chronology of abandoned loops at each reverse migration time step

At a given reverse migration time step, the paleo-channel is susceptible to capture n abandoned meanders, which must be selected among the most recent members of each group of abandoned loops (Fig. 7). For this, we compute probabilities of integration for the candidate abandoned meanders given the current channel path and based on geometric considerations.

Choosing the number n may be achieved by studying natural, flume tank or numerical analogs. The Appendix A presents such an example, where the Mississippi database is used to assess the number of abandoned meanders per unit length and per year.

A number of abandoned meanders is sampled from the Poisson probability distribution of cutoff rate. If the sampled number is non-null, all abandoned meanders with a relative age equal to 0 are candidates for integration at the current time step. They are as numerous as the number of observed groups of abandoned meanders. A likelihood of integration is then computed for each considered abandoned loop.

Once the number of abandoned meanders to integrate at the current time step has been drawn, the likelihood of integration of each candidate abandoned meander needs to be evaluated. In the forward sense, the process of abandonment by stream connection between two bends of a meander loop implies that the loop is close and roughly orthogonal to the main channel trajectory at the time of abandonment (e.g., [START_REF] Lonsdale | Cut-offs at an abyssal meander south of Iceland[END_REF][START_REF] Erskine | Alluvial cutoffs as indicators of former channel conditions[END_REF]. Indeed, the cutoff can be due to the proximity of the two bends (neck cutoff) or to the development of a shortcut after the modification of the hydrodynamic flow (chute cutoff) [START_REF] Gay | Evolution of cutoffs across meander necks in Powder River, Montana, USA[END_REF][START_REF] Constantine | A mechanism of chute cutoff along large meandering rivers with uniform floodplain topography[END_REF][START_REF] Constantine | Sediment supply as a driver of river meandering and floodplain evolution in the Amazon Basin[END_REF]. The latter is more common in fluvial settings than in turbidite ones. In both cases, the meander loop is by-passed by the modification of the channel flow through the neck or chute cutoff of the meander loop. In the proposed reverse migration method, we simply use proximity and The likelihood is influenced by the orientation and proximity of the paleomeander relative to the main channel path. In terms of orientation, an abandoned meander that is perpendicular to the main channel has the highest chance to be integrated. Conversely, a meander parallel to the main channel path is not accepted in our approach because abandonment is considered as improbable and because its geometric integration of the meander would be delicate [START_REF] Parquer | Reconstruction of channelized systems through a conditioned reverse migration method[END_REF].

• The nearest distance is defined as the shortest path between the segment apex-cutoff and the main channel path (Fig. 8). An abandoned meander located too close to or too far from the main channel path has no chance to be considered as abandoned by the current channel. For integration to be possible, the domain is therefore limited so that the distance is less than half the wavelength of the channel, but greater than one-eighth of the half-wavelength of the channel. This value has been chosen to avoid problematic proximity between the abandoned meander and the channel path. Within these bounds, a linear decrease dependent on the channel half-wavelength HW represents the likeli-

O i P distance (i) = 0 if d i < HW 8 or d i > HW 1 - d i -HW 8 7HW 8 if HW 8 ≤ d i ≤ HW , (1)
in which d i is the shortest distance in meters between the abandoned meander O i and the main channel path.

• The abandoned meander can be oriented by its two extremities and its apex (Fig. 9). The angle α i formed by this direction and the main channel direction closest to the tips of the abandoned meander is used to define the likelihood of integration

P orientation (i) = 0 if α i < 15 1 -90-α i 75 if α i ≥ 15 , (2)
By default, the half-wavelength HW is estimated from the current channel path. It could also be deduced from an analog [START_REF] Inglis | The behavior and control of rivers and canals[END_REF][START_REF] Leopold | River meanders[END_REF][START_REF] Williams | River meanders and channel size[END_REF][START_REF] Flood | Quantitative characteristics of sinuous distributary channels on the Amazon Deep-Sea Fan[END_REF][START_REF] Malavoi | River hydrogeomorphology: A primer, Onema edn Miall AD[END_REF]. Finally, both likelihoods are combined to rank the abandoned meanders according to their integration likelihood:

P integration (i) = P distance (i) + P orientation (i) 2 . ( 3
)
The above likelihoods are then transformed into a probability of integration of the abandoned loop O i by normalizing over all the candidate loops O j of minimal age:

P integrationnormalized (i) = P integration (i) P integration (j) . ( 4
)
see an example on in Fig. 10. Then, the n paleo-meanders to integrate at this time step are randomly sampled from this probability distribution (Fig. 10).

Application to paleo-channel reconstruction

The technique described above makes it possible to find among the observed abandoned meanders, at a given time step, the best candidates for integration into the path of the current channel. However, migration is a continuous process and after the integration of a first group of paleo-meanders, the main channel path has to be updated before considering the remaining meanders. In this paper, we propose to achieve this paleo-reconstruction with the ChaRMigS reverse migration method proposed by [START_REF] Parquer | Reconstruction of channelized systems through a conditioned reverse migration method[END_REF]. This reverse migration technique permits us to simulate the paleotrajectories of a given channelized system starting from a particular channel path interpreted from available data, in addition to probability distributions of migration amplitudes and abandoned meander paths. As illustrated by Fig. 3, enriching the reverse channel migration with meander integration and the local relative age update enriches the production of multiple possible global chronologies and associated 3D system architectures.

Integration in the reverse migration simulation

When a paleo-meander has been selected, its integration is performed at the next reverse migration time step (Fig. 11). Once the paleo-meander has been integrated into the main channel path, its integration time step is registered in the global chronology, and it is removed from the integration waiting list. The older abandoned loops remain in the corresponding group and their local chronologies are decremented, so that they can be considered in the following iterations. For instance, in the channel belt shown in Fig. 11, when the entity A is selected for abandonment, the relative ages for meander loops B and C are set to 0 and 1, respectively.

In the above workflow, the integration of an abandoned meander loop in the main channel path may fail [START_REF] Parquer | Reconstruction of channelized systems through a conditioned reverse migration method[END_REF], for example, if the abandoned meander tips make a large angle with the current channel, or when an overlap with the main channel path can be observed. In this case, the integration of the abandoned loop is postponed to the next reverse migration step.

Ideally, the best would be to let the algorithm run until all paleo-loops are chronologically ordered and integrated into the main channel path, but integration failure and postponement may lead to impractical simulation times. Therefore, this workflow is repeated until the maximal number of time steps chosen for the realization is reached (Fig. 12).

The simulation of the global chronology of paleo-meanders is closely related to the reverse migration of the complete system. Both simulations are run simultaneously with an update of all distance and orientation computations at each time step so that probabilities of integration are consistent with the current channel path. This facilitates the paleo-meander integration while respecting the local chronologies. However, the strict enforcement of local chronologies may have an adverse effect in terms of honoring all observations. Indeed, if the younger abandoned meander of a group fails to be integrated, the other abandoned meanders of the same group will never be considered for integration. This suggests that the reverse migration method may generate channel paths which are not compatible with some abandoned meander interpretations. As reverse migration (without abandoned meander integration) essentially tends to decrease channel curvature, this may occur because some abandoned meanders are missing from the interpretation, or because the channel path used to initiate the algorithm is wrong. Running several simulations may help assessing whether such inconsistencies appear systematically or only in some realizations.

McMurray Formation application

We applied the proposed method, to a 3D seismic stratal slice image of the McMurray Formation in Alberta (Canada) (Fig. 13a). Locally, the Mc-Murray Formation occurs at a depth of 350 to 500 m below the surface, and here the 3D seismic image acquisition covers almost 600 km 2 (Fig. 13a). The Lower Cretaceous McMurray Formation in northeastern Alberta (Canada) is a record of a complex continental to estuarine depositional history [START_REF] Mossop | Deep channel sedimentation in the Lower Cretaceous McMurray Formation, Athabasca oil sands, Alberta[END_REF][START_REF] Ranger | Elements of a stratigraphic framework for the McMurray Formation in south Athabasca area, Alberta[END_REF][START_REF] Hein | The Athabasca oil sands a regional geological perspective, Fort McMurray area, Alberta, Canada[END_REF]. On the basis of studies on 3D seismic data, scaling relationships, and detritalzircon geochronology, the hypothesis that the catchment for these fluvial sediments was of continental scale, similar to that of the modern Mississippi River has been made [START_REF] Musial | Subsurface and outcrop characterization of large tidally influenced point bars of the Cretaceous McMurray Formation (Alberta, Canada)[END_REF]Benyon et al., 2014;[START_REF] Blum | Mid-Cretaceous to Paleocene North American drainage reorganization from detrital zircons[END_REF][START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF]. Evidence of tidal influence has been deduced in the fluvial point bar-dominated deposits within the study area (e.g., [START_REF] Hubbard | Seismic geomorphology and sedimentology of a tidally influenced river deposit, Lower Cretaceous Athabasca oil sands, Alberta, Canada[END_REF][START_REF] Labrecque | Sedimentology and stratigraphic architecture of a point bar deposit, Lower Cretaceous McMurray Formation, Alberta, Canada[END_REF][START_REF] Musial | Subsurface and outcrop characterization of large tidally influenced point bars of the Cretaceous McMurray Formation (Alberta, Canada)[END_REF]. Thus, considering seismic data, zircon geochronology, and scaling laws, both Mississippi and McMurray systems are interpreted to be quite similar in terms of environmental context and system dimensions [START_REF] Durkin | Threedimensional reconstruction of meander-belt evolution, Cretaceous McMurray Formation, Alberta Foreland Basin, Canada[END_REF]. Therefore, we extracted the frequency of abandonment from the Mississippi analog (Appendix A and Table 2) here approximated by a Gaussian probability distribution. The cutoff rate has been computed in years from the Mississippi database. The rate is transformed for its use in the McMurray Formation in terms of time steps. As the conversion of years in time step is not obvious, in this study, the cutoff rate for each time step has been taken as the computed cutoff rate of the most recent survey (Appendix A and Table 2). In the same way, migration offset probability distributions have been inferred by the Mississippi analog [START_REF] Parquer | Reconstruction of channelized systems through a conditioned reverse migration method[END_REF].

The choice of this seismic image has been motivated by the high number of abandonment shapes that were observable and by the numerous studies made on this formation, making comparisons possible. However, the resolution level of this image is significantly better than the one targeted by the designed method. The small details such as scroll bars are not considered here: thus, in this study, we only used the ten clear abandoned meander paths that were observed and digitized (Fig. 13b). Nine intersections were brought out by the semi-automatic analysis of the channel belt (Fig. 13b). These intersections divide the abandoned meanders into four distinct groups (Fig. 13c). Then, inside each of the four groups, the cross-cutting relationships between paleo-meanders and the associated interpretation of relative ages led us to determine local chronologies (Fig. 13d).

Thanks to the stochastic reverse migration results (Figs. 14b1 and 14b2 and 14b3 and 14b4), we found various global chronologies (Figs. 14c1 and 14c2 and 14c3 and 14c4). All of them respect the observed truncation erosive patterns and the interpreted relative chronologies. The variability brought by the stochasticity of the reverse migration favors differently the integration of the various abandoned meanders through their abandonment probability, which is linked to the updated channel path.

We can compare this work with that of [START_REF] Durkin | Threedimensional reconstruction of meander-belt evolution, Cretaceous McMurray Formation, Alberta Foreland Basin, Canada[END_REF] in which one reconstruction has been proposed for the complete system. The reconstruction of Durkin ( 2016) is based on precise analyses of the observed point bars in order to chronologically classify and link them relatively one to another. This leads to a precise quantification of preservation rates deduced for each area and time lapses [START_REF] Durkin | Evolution of fluvial meander-belt deposits and implications for the completeness of the stratigraphic record[END_REF]. In contrast, our simulation results permit us to span the range of possible architectures from less information, as we only used sparse meander loop picks and much less detailed cross-cutting analyses of structures on the seismic images. Because we never access to the absolute ages of the structures, such a detailed study on the preservation rate like the one made by [START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF] is not possible in our approach. However, even though less information was used, we can observe that the relative ages obtained on the four proposed realizations are almost the same as those deduced by [START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF] (Table 1). Among the differences, for example, the abandoned meander C located on the east bank of the channel, according to [START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF], should be integrated third, after the abandoned meanders J and G. In our simulation, none of the realizations produced the exact same output. This could be explained by its shape, which is closer to a chute cutoff than to a neck cutoff, and for which the probabilistic criteria of proximity and orientation with respect to the channel are therefore less relevant. Nevertheless, it would be interesting to run other simulations in order to see if this scenario is rare. Another difference remains in the abandoned meander A which, should, according to the extensive study of [START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF], be among the first ones to be integrated. In our approach, the closest and the most perpendicular one of the meander belt is expected to have the highest probability of integration. This illustrates that the proximity to the channel is not necessarily the only criterion to be used for deciding about abandoned meander integration. For example, if available, depth should also be considered. A study could be done by imposing the chronology order and comparing the reverse migrated paths to the ones reconstructed by [START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF].

Another interesting point is to compare the sinuosity variations of [START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF] and the simulation results of our method. [START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF] proposed a variation of the sinuosity between 1.11 for the last observed (youngest) channel path and 1.68 for the oldest reconstructed one with two peaks at 2.61 and 2.52. The four realizations presented in Fig. 14 present a sinuosity between 1.33 for the last (youngest) channel path and between 1.20 and 2.08 for the oldest reconstructed one with peaks between 2.35 and 3.12 depending of the considered realization. The average sinuosity in the reconstruction of [START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF] is equal to 2.08 versus 1.45 for the four proposed stochastic reconstructions. Overall, these sinuosity figures are comparable for both methods and are in the range of sinuosities observed in the modern Mississippi analog.

Discussion

In the proposed method, model consistency is driven by the geometric evolution in reverse migration steps, under strict relative age constraints obtained by cross-cutting analysis or other interpretation criteria. The gradual evolution is key (e.g., [START_REF] Sylvester | A model of submarine channellevee evolution based on channel trajectories: Implications for stratigraphic architecture[END_REF][START_REF] Parquer | Reconstruction of channelized systems through a conditioned reverse migration method[END_REF] to gain global consistency. Thus, integrating this chronology technique inside the reverse migration algorithm plays a major part in moving toward a better system reconstruction. Multiple possible global chronologies can be obtained by the proposed method. All of them respect local chronologies and permit the production of a range of possible geometries of architectural elements through stochastic samplings at diverse steps of the workflow. The different sets of possible chronologies translate the inherent uncertainties due to the incompleteness of the sedimentary record, which is often exacerbated by incomplete or fuzzy observations. In this paper, we focus on a scenario where the abandoned meanders and the most recent channel path are interpreted on satellite, seismic, or LIDAR images. In its current implementation, the method relies on the expert-based extrapolation of eroded or non-visible parts of abandoned meanders before reverse migration of the complete system (Fig. 1 andAppendix A). Namely, the method requires the complete reconstruction of the abandoned meander from one extremity to the other even when it is not observable. Ideally, it would be good to avoid such a manual extrapolation, which may introduce a bias in the results. Indeed, the extrapolated extremities may impact the relative ages obtained by intersection analyses. Also, these manual extrapolation of eroded parts certainly affect the probabilities of integration and whether the actual meander integration is successful or not. For the sake of the method's efficiency, we can imagine that extremities would be simulated just before the integration stage so as to enhance the probability of success of the integration. Our workflow could also use interpretations generated or simulated by computational methods (e.g., image analysis (e.g., [START_REF] Ruiu | Semi automatic interpretation of 3D sedimentological structures on geologic images: An object-based approach[END_REF] or simulation of unseen or eroded channel parts (e.g., [START_REF] Pyrcz | ALLUVSIM: A program for eventbased stochastic modeling of fluvial depositional systems[END_REF][START_REF] Rongier | Stochastic simulation of channelized sedimentary bodies using a constrained L-system[END_REF]). The simulation of extremities for all candidate abandoned meanders would add another source of stochasticity to the algorithm. Also, this simulation could be done differently if the abandoned meander is supposed to be a neck cutoff or a chute cutoff. This could help to address the chronological differences of meanders C and J discussed in Section 4.2 in the proposed method. The abandonment ages obtained after the final reconstruction step are uncertain, both because of uncertainties about migration rates (which affect the actual time step duration) and because the stochastic method can produce several ages for the same abandoned meander depending on the realization. This translates the ability of the method to sample chronological uncertainties about meander abandonment. These aspects could be further studied by propagating time step and migration rate uncertainties during the process to generate a more robust probability distribution of abandonment ages.

In this workflow, we also assume that all abandoned meanders have been identified. However, some abandoned meanders can be completely eroded, as suggested by the cutoff frequency study made on the Mississippi (Ap-pendix A). Some others may not be observable due to low imaging resolution. This may be one reason why applications of our reconstruction method to real data could differ from results obtained by other techniques. Indeed, in contrast to our approach, the manual reconstructions have the freedom to interpret an abandoned loop for which no trace remains in order to fit the interpreter's observations. In our methodology, the reverse migration without meander integration tends to decrease channel curvature and may yield isolated abandoned meanders that are too far from the channel to be integrated. Considering fully eroded or below resolution abandoned meanders could certainly help to avoid this problem and facilitate the integration of all the visible meanders. Moreover, the cutoff rate frequencies inferred from the Mississippi data suggests that abandoned meanders are no longer observable because of erosion. Thus, additional abandoned meanders could be simulated to reproduce more realistic past channel evolution. Such analog study highlights the need of abandoned meander simulation in reconstruction purposes of channelized systems. In the case of recurrent non-integrated loops, it could be interesting to draw or simulate missing abandoned loops in some particular locations to test the impact of meander simulation on the integration rate.

The probability distribution of the spatio-temporal rate of meander abandonment is crucial for the reconstitution of the complete system deposits. Indeed, this probability impacts the number of reverse time steps between the integration of two abandoned meanders and thus the architecture and the connectivity of the simulated channel paths. Such information is almost impossible to extract from seismic or satellite images. A general solution to choose appropriate reverse migration and meander abandonment rate parameters is to use analogs (e.g., outcrops, modern systems or processbased models) [START_REF] Schwenk | Meander cutoffs non locally accelerate upstream and downstream migration and channel widening[END_REF][START_REF] Sylvester | High curvatures drive river meandering[END_REF]. These parameters are also controlled by external factors such as water and sediment discharge [START_REF] Constantine | A mechanism of chute cutoff along large meandering rivers with uniform floodplain topography[END_REF], topography, tectonic settings, subsidence or climate. For example, [START_REF] Durkin | Threedimensional reconstruction of meander-belt evolution, Cretaceous McMurray Formation, Alberta Foreland Basin, Canada[END_REF] observed on his reconstruction of the McMurray Formation that the migration was in a first period mainly lateral with expansion of meanders and in a second phase more important in the downstream direction. The explanation mainly relies on the hypothesis that the shale infills of abandoned meanders create a narrower space for migration of the main channel path as time goes by. Such observation could be integrated in the proposed method by tuning the relative importance of both horizontal migration direction time step by time step [START_REF] Parquer | Reconstruction of channelized systems through a conditioned reverse migration method[END_REF]. More generally, rules on the evolution of influencing factors could be translated by varying the parameters of the various probability distributions used in the reverse migration.

As the complete method ultimately reconstructs possible temporal evolution time steps of the channel from available evidence, it would also be useful to explicitly model the geometric layout of depositional heterogeneities based on paleo-channels [START_REF] Yan | A 3D forward stratigraphic model of fluvial meander-bend evolution for prediction of pointbar lithofacies architecture[END_REF]. The consideration of heterogeneities is of first importance for flow simulation in order to better understand the reservoir architecture [START_REF] Jacquemyn | Surface-based geological reservoir modeling using grid-free NURBS curves and surfaces[END_REF]. However, these heterogeneities are of highly variable size and extension: low angle contacts and relatively small features between architectural elements could lead to meshing issues prior to flow simulation [START_REF] Mustapha | An efficient method for discretizing 3D fractured media for subsurface flow and transport simulations[END_REF][START_REF] Pellerin | Automatic surface remeshing of 3D structural models at specified resolution: A method based on Voronoi diagrams[END_REF].

Our approach is applicable to most settings where image quality is relatively low. The method merely needs the digitizations of the main channel path and remaining loops as well as various probability distributions. The scroll bar patterns, only observable in high precision images, are not considered. Honoring them is, indeed, mainly linked to the correct inference of the reverse migration parameters, which is beyond the scope of this paper. One of the main drawbacks is that the reconstructed paleo-paths do not always correspond to the actual paleo-paths. This is precisely the goal of stochastic simulation to approach the truth by a range of possible realizations of reconstructed geometries. However, when some internal structures are observable such as point bars on some areas, we could reject the realizations when the angle between the reconstructed paleo-trajectory deviates from the scroll bars. More immediately, relative ages of abandoned meanders deduced from scroll bar truncations could serve as additional constraints of our method. Another improvement could be made by repulsing the channel reverse migration away from areas where an abandoned meander has been observed. For the moment, a simple rejection of the realizations presenting such an unlikely result is possible. These improvements could be welcomed in the future, but could bear a significant impact on the computational performance of the method. At present, the speed and memory usage of the proposed method is already relatively large. Indeed, the update of the main channel path at each time step in order to compute updated probabilities of integration of each abandoned meander has a direct impact on the algorithm's efficiency. Thus, generating a realization such as the ones presented in Fig. 14 takes about two hours on a 64-bit system with a 2.00 GHz processor Intel Core i7-4510U and 16 GB of RAM. Optimizations of our implementation are probably possible, but parallelization is made difficult by the sequential nature of the reverse migration.

More fundamentally, one may question the value of an automatic method to reach the level of interpretation of an expert. We believe that the promise of automation is mainly to generate several possible interpretations based on expert input. This may not be very useful when images are very accurate and allow an interpreter to directly see all the features of interest and to come up with a consistent evolution. However, it can be critical in most cases where ambiguity exists in the data. In these cases, it has been clearly shown that interpretation bias exists in structural studies [START_REF] Bond | What do you think this is? "conceptual uncertainty" in geoscience interpretation[END_REF][START_REF] Bond | Uncertainty in structural interpretation: Lessons to be learnt[END_REF]. We are not aware of similar studies on bias in channel interpretation, but it is likely to also exist for psychological reasons [START_REF] Kahneman | Thinking, fast and slow[END_REF]. In any case, coming up with several possible scenarios that can be scrutinized by experts is certainly an interesting way forward to understand subsurface uncertainty.

Conclusion

We proposed a stochastic semi-automated tool to determine the possible relative ages of abandoned meander loops and reconstruct the past geometric evolution of the associated channels. This approach is useful when a precise study of the channel belt (e.g., paleo current, point bars) or absolute dating method is not possible. Indeed, this tool overcomes these difficulties by considering relative age constraints between abandoned meanders and by simulating stochastically the global chronology by a reverse migration process and geometric rules. It could also be used as an assistance for more detailed, expert-based work such as [START_REF] Durkin | Threedimensional reconstruction of meander-belt evolution, Cretaceous McMurray Formation, Alberta Foreland Basin, Canada[END_REF] or [START_REF] Hubbard | Seismic geomorphology and sedimentology of a tidally influenced river deposit, Lower Cretaceous Athabasca oil sands, Alberta, Canada[END_REF], by proposing quick stochastic realizations of possible 3D geometries of a paleosystem.

The method provides possible global ages for the abandoned loops based on geometric considerations. This age significantly helps the reverse migration method finding alternative channel evolution paths as compared to the reverse migration method of [START_REF] Parquer | Reconstruction of channelized systems through a conditioned reverse migration method[END_REF]. It also increases the reliability of the reconstructed geometries by identifying at each time step the most likely abandoned meanders to be integrated between reverse migration steps. The application on the McMurray seismic image has shown that, with fewer observations, the relative meander ages we simulated are quite consistent with those obtained by [START_REF] Durkin | Threedimensional reconstruction of meander-belt evolution, Cretaceous McMurray Formation, Alberta Foreland Basin, Canada[END_REF] (Table 1). The sinuosity of the system is also comparable. Applications in better known systems such as the Mississippi could help assess the ability of the method to recover the correct paleo-meander trajectories. However, the exact ages of the reconstructed channels will only be probabilistic, as long as the migration rates cannot be exactly constrained and the global uncertainty impossible to wipe out.

One limit of this technique is that the result is only conditioned by the observed paleo-channels. However, many of the abandoned meander loops of the system did not survive later evolution of the system. Simulating eroded abandoned loops is an ongoing work which could help to better extend the reverse migrated system to the complete meander belt and, so, to assess the potential geometries of the reservoir and the associated heterogeneities.

We would like to apply our method to other seismic data from turbidite or fluvial contexts. Using 3D seismic data would be interesting to test the conditioning to paleo-meanders at various depths. Differences in fuzziness depending on the channel base depth due to sub-seismic level of detail could also be a subject of investigation.

Long-term future works concern going toward validation of the proposed model geometries through the simulation of the heterogeneities inside the channelized system and the calibration to hydrogeological or reservoir production data. many kinds. However, the Mississippi cannot be considered as a reliable analog for any kind of channelized system, such as for example a highly tidally influenced system or a strict marine system.The last four paths of the Mississippi river mapped with precision by [START_REF] Fisk | Geological investigation of the alluvial valley of the lower Mississippi River[END_REF] have been digitized along 1200 kilometers (Fig. 15). They span more than 250 years of river evolution.

The mapping of the cutoffs is more or less clear. Indeed, sometimes, the flow continues to pass in the meander loop that have been shortcut. Such a configuration is illustrated by the southern cutoff on the last channel path in Fig. 15. This is translated in the map by a divergence of the current flow at the location of the cutoff.

On these Mississippi paleo-channel digitizations, we counted the clear cutoffs (without divergent contemporary flows) between each pair of successive paths. We considered the divergent branching cutoff in the standard deviations of each period. The results are presented in Table 2.

The differences in the number of observed abandoned meanders created between pairs of digitized paths can be explained by the decrease over time of the preservation probability of structures (e.g., [START_REF] Sadler | Sediment accumulation rates and the completeness of stratigraphic sections[END_REF]. This study gives us some minimal rates assuming that the driving parameters remained the same. In order to estimate the real number of abandoned meanders that are created per year in an analogue such as for example in the McMurray Formation in Section 4.2, we should take the most recent survey, so a mean

Figure 1 :

 1 Figure 1: Constraints detected from image analysis depending on the image resolution: channel and abandoned meander paths, cross-cutting patterns and point bars (a) LIDAR data of the Siiponjoki, Finland, [Geological survey of Finland 64 • 12 33.28"N 23 • 45 31.56"E] (b) Satellite image of the Wabasca river, Canada, [Google Earth 58 • 11 10.59"N 115 • 24 42.00"O] (c1) High resolution seismic image of McMurray Formation, Alberta, Canada after Durkin et al. (2017), from ConocoPhillips Canada. (c2) The analysis of internal point bar structures permits to analyze the gradual evolution of the channel. The relative ages of abandoned meanders can only be determined when cross-cutting by other paleo-channel paths is observed.

Figure 2 :

 2 Figure 2: Global workflow to deduce the 3D channelized system architecture and the global chronology of abandoned meanders from geometric channel belt analysis

Figure 3 :Figure 4 :

 34 Figure 3: Workflow of the Channel Reverse Migration Simulation (ChaRMigS). The integration of the chronology simulation tool is highlighted in red. Instead of an abandonment age equal to a fixed value or to a value sampled from some a priori age probability distribution Parquer et al. (2017), this paper proposes a simulation of the chronology for each time step of the reverse migration (Algorithm 1), which accounts for the current channel geometry. Modified from[START_REF] Parquer | Reconstruction of channelized systems through a conditioned reverse migration method[END_REF]

Figure 5 :

 5 Figure 5: (a) Automatic detection of intersections between interpreted paleo-meanders. (b) Detection of groups of abandoned loops through automatic detection of intersections

Figure 6 :

 6 Figure 6: Deduction of local chronologies of meanders labeled A to G inside groups through relative ages observation. The abandoned meanders not eroded by another abandoned meander are automatically candidate for abandonment at the current time step. The eroded ones have an added value of local chronology and necessarily correspond to an anterior path of the main channel. The four yellow abandoned meanders are all first youngest in their respective group but cannot be ordered only by the erosive relationships between abandoned meanders.

Figure 7 :

 7 Figure 7: Workflow of abandoned meander draw at a given reverse migration time step

Figure 8 :

 8 Figure 8: Influence of the abandoned meander distance inside its meander belt on its integration probability P distance , with HW the half-wavelength of the channel path (a) likelihood of integration according to the abandoned meander distance (b) computation of the shortest distance of the abandoned meander and the main channel path

Figure 9 :

 9 Figure9: Influence of the abandoned meander orientation inside its meander belt on its integration probability P orientation (a) likelihood of integration according to the abandoned meander orientation (b) computation of the orientation of the abandoned meander and the main channel path

Figure 10 :

 10 Figure10: The abandoned meanders are drawn among the youngest ones (in dark blue in the table, in light yellow in the schema, the ones with the relative chronology value equal to zero) according to their probability of integration (in red)

Figure 11 :

 11 Figure11: Update of the lists of the abandoned meanders which are still candidate for abandonment after one of them has been drawn at the previous time step. The entity B has now a local chronology value equal to zero. However, its abandonment probability is also equal to zero because its location is too close to the updated channel path. The same configuration is observed for the abandoned meander G.

Figure 12 :

 12 Figure12: The reverse migration conditioned to paleo-meanders proceeds until the number of time steps is reached. At the end of the process, as many as possible paleo-meanders have been integrated according to their abandonment probabilities and a global chronology is deduced.

Figure 13 :

 13 Figure 13: Application to the McMurray Formation, Alberta, Canada. (a) Seismic stratal slice image after Durkin et al. (2017), from ConocoPhillips Canada. (b) Digitization of the ten observed abandoned loops (labeled A to J) and automatic detection of intersections between digitized paths (highlighted by red crosses). (c) Automatic detection of groups through intersections between digitizations. (d) Local chronologically ordering inside groups

Figure 14 :

 14 Figure 14: Application to the McMurray Formation, Alberta, Canada. (a) Semi-automatic chronological ordering inside groups. (b1), (b2), (b3) and (b4) Four realizations of reverse migration (100 migration time steps, probability distribution of amplitude horizontal offset N(1000, 300) in meters). (c1), (c2), (c3) and (c4) Four resulting realizations of global chronology simulation processed during the reverse migration

Figure 15 :

 15 Figure 15: Picking of the four last trajectories of the Mississippi in 1944 (in white), 1880 (in green), 1820 (in red) and 1765 (in blue) on a map [http://lmvmapping.erdc.usace.army.mil/]

migration step IF abandonment age is reached definition of the migration directions FOR EACH half- meander Z axis

	INPUT
		FOR EACH
			reverse
	-N number of
	wanted reverse
	migration steps
	-H & Z migration
	offset distributions
	-abandoned loop digitization &
	relative abandonment age estimation
	IF ANY
	ABANDONED LOOP INTEGRATION	REVERSE MIGRATION	OUTPUT -channel belt geometry due to channel stacking

Table 1 :

 1 Summary of the global relative chronologies found by the study of[START_REF] Durkin | The evolution of fluvial meander belts and their product in the rock record[END_REF] and our four realizations.

now at Geological Survey of Canada, Ottawa, Canada

Acknowledgments

This work was performed in the frame of the RING project (http://ring.georessources.univlorraine.fr/) at Université de Lorraine. We would like to thank for their support the industrial and academic sponsors of the RING-GOCAD Consortium managed by ASGA. The authors would also like to thank Paul Durkin for sharing the McMurray Formation seismic image (courtesy of ConocoPhilips Canada). Software corresponding to this paper is available in the GoNURBS plugin of SKUA-Gocad. We also acknowledge Paradigm for providing the SKUA-Gocad Software and API. The authors thank Paul Durkin and an anonymous rewiever for their constructive remarks that helped us to improve this paper.

Benyon C, Leier A, Leckie DA, Webb A, Hubbard SM, Gehrels G (2014) Provenance of the Cretaceous Athabasca oil sands, Canada: Implications for continental-scale sediment transport. Journal of Sedimentary Research 84(2): 136-143, doi: 10.2110/jsr.2014.16

Appendices A. Cutoff rate determination based on Mississippi analog study

The richness and the public availability of such a big system makes the Mississippi an interesting analog for extracting probability distributions of