Multiple-point Statistics Simulation Research Articles

Generative Adversarial Networks to incorporate the Training Image uncertainty in multiple-point statistics simulation

Multiple-point geostatistics (MPS) algorithms have been extensively used for geostatistical reservoir modeling. These techniques are known for building geostatistical simulations with curvilinear features, which are crucial for applications in flow simulation. One challenge in using MPS algorithms is obtaining a representative Training Image (TI). The result is that often one TI is used to condition the entire set of geostatistical simulations, and the TI uncertainty is disregarded. If the TI is uncertain, this uncertainty should be considered in the geostatistical simulations. Otherwise, the spatial uncertainty represented by the geostatistical simulations is underestimated. Incorporating the TI uncertainty requires building a data set or catalog of plausible TIs. In this context, we propose to use Generative Adversarial Networks (GANs) to build this catalog of TIs prior to the MPS simulation. Each TI generated by GAN is used to condition one realization built by the chosen MPS algorithm. The proposed workflow combines the strength of GAN to generate TIs with the strength of MPS algorithms to consider several sources of information to condition the simulations. A case study is shown, where the methodology is compared against the traditional MPS workflow, which uses a single TI. The results show that the proposed method successfully incorporated the TI uncertainty into the geostatistical workflow. The proposed method built numerical models with higher uncertainty and variability and, thus, provides a more realistic quantification of the space of uncertainty.

Markov chain Monte Carlo for seismic facies classification

Seismic facies classification aims to predict a facies model, or a set of facies models, from measured seismic data. We focus on stochastic classification methods to estimate the probability distribution of facies conditioned on seismic data. Bayesian classification methods based on analytical solutions are generally applied to seismically inverted elastic attributes or petrophysical properties rather than measured seismic data, and they do not include spatial correlation models in the prior distribution. Instead, iterative stochastic methods, such as Monte Carlo acceptance-rejection sampling and approximate Bayesian computation, can be directly applied to seismic data and can be implemented with complex prior models; however, their convergence is generally slow because it requires the simulation of a large number of facies models from the prior distribution. We develop an efficient implementation of a Markov chain Monte Carlo (MCMC) approach that adopts complex prior models, such as multiple-point statistics simulations based on a training image, to generate geologically realistic facies realizations. The novelty of the approach is that the proposal distribution of the proposed MCMC method is based on a gradual deformation algorithm, in which the proposal distribution is expressed as a linear combination of the indicator variables of the previously accepted model and a random prior simulation, according to a uniform random weight. The synthetic examples indicate accurate facies predictions with a success rate of approximately 80% for well-log facies classification based on elastic attributes and approximately 60% for seismic facies classification based on seismic data. The inversion is compared to an MCMC method with prior models sampled from a first-order Markov chain and Bayesian facies classification.

Reconstructing Three-dimensional geological structures by the Multiple-point statistics method coupled with a deep neural network: A case study of a metro station in Guangzhou, China

A credible large-scale three-dimensional (3D) geological model constructed from a small amount of local high-density data provides a vital data infrastructure for the design and construction of metro lines. The multiple-point statistics (MPS) method has successfully generated 3D geological models in many fields. However, the global cognition of the spatial relations between geological objects is difficult to reconstruct using a local optimization strategy with moving templates in the MPS simulation. This study presented an MPS simulation approach for constructing 3D geological structures with two-dimensional (2D) cross-sections, where the MPS process was coupled with a fully connected deep neural network (FCDNN) simulation. The core idea of the proposed method was that the global features of each geological object were generated from the trained FCDNN, the input parameter of which was the subsurface depth. Then, the entire model was constructed on the initial model using the following iterative MPS process: the semantic relationship of geological objects was extracted from 2D cross-sections and used as a constraint in the modeling process, which ensured the order of the stratigraphic sequence and rationality of the geological structure in the final model. The joint objective function of the proposed method was composed of the loss function of the FCDNN and conventional MPS objective function constrained by the spatial semantics of geological objects. The concrete example in the Jiangtai Road Metro Station of Guangzhou Line 11, China illustrated that the details of the complex geological structures could be obtained accurately by the presented method. The global features from the trained FCDNN provided the basic geometry for the final model, which overcame the limitations of the conventional MPS method.

Hydrogeological multiple-point statistics inversion by adaptive sequential Monte Carlo

For strongly non-linear and high-dimensional inverse problems, Markov chain Monte Carlo (MCMC) methods may fail to properly explore the posterior probability density function (PDF) given a realistic computational budget and are generally poorly amenable to parallelization. Particle methods approximate the posterior PDF using the states and weights of a population of evolving particles and they are very well suited to parallelization. We focus on adaptive sequential Monte Carlo (ASMC), an extension of annealed importance sampling (AIS). In AIS and ASMC, importance sampling is performed over a sequence of intermediate distributions, known as power posteriors, linking the prior to the posterior PDF. The AIS and ASMC algorithms also provide estimates of the evidence (marginal likelihood) as needed for Bayesian model selection, at basically no additional cost. ASMC performs better than AIS as it adaptively tunes the tempering schedule and performs resampling of particles when the variance of the particle weights becomes too large. We consider a challenging synthetic groundwater transport inverse problem with a categorical channelized 2D hydraulic conductivity field defined such that the posterior facies distribution includes two distinct modes. The model proposals are obtained by iteratively re-simulating a fraction of the current model using conditional multiple-point statistics (MPS) simulations. We examine how ASMC explores the posterior PDF and compare with results obtained with parallel tempering (PT), a state-of-the-art MCMC inversion approach that runs multiple interacting chains targeting different power posteriors. For a similar computational budget, ASMC outperforms PT as the ASMC-derived models fit the data better and recover the reference likelihood. Moreover, we show that ASMC partly retrieves both posterior modes, while none of them is recovered by PT. Lastly, we demonstrate how the power posteriors obtained by ASMC can be used to assess the influence of the assumed data errors on the posterior means and variances, as well as on the evidence. We suggest that ASMC can advantageously replace MCMC for solving many challenging inverse problems arising in the field of water resources.

Ice volume and basal topography estimation using geostatistical methods and ground-penetrating radar measurements: application to the Tsanfleuron and Scex Rouge glaciers, Swiss Alps

Abstract. Ground-penetrating radar (GPR) is widely used for determining mountain glacier thickness. However, this method provides thickness data only along the acquisition lines, and therefore interpolation has to be made between them. Depending on the interpolation strategy, calculated ice volumes can differ and can lack an accurate error estimation. Furthermore, glacial basal topography is often characterized by complex geomorphological features, which can be hard to reproduce using classical interpolation methods, especially when the field data are sparse or when the morphological features are too complex. This study investigates the applicability of multiple-point statistics (MPS) simulations to interpolate glacier bedrock topography using GPR measurements. In 2018, a dense GPR data set was acquired on the Tsanfleuron Glacier (Switzerland). These data were used as the source for a bedrock interpolation. The results obtained with the direct-sampling MPS method are compared against those obtained with kriging and sequential Gaussian simulations (SGSs) on both a synthetic data set – with known reference volume and bedrock topography – and the real data underlying the Tsanfleuron Glacier. Using the MPS modeled bedrock, the ice volume for the Scex Rouge and Tsanfleuron glaciers is estimated to be 113.9 ± 1.6 million cubic meters. The direct-sampling approach, unlike the SGS and kriging, allowed not only an accurate volume estimation but also the generation of a set of realistic bedrock simulations. The complex karstic geomorphological features are reproduced and can be used to significantly improve for example the precision of subglacial flow estimation.

Hierarchical MPS-Based Three-Dimensional Geological Structure Reconstruction with Two-Dimensional Image(s)

Multiple-point statistics (MPS) is a useful approach to reconstruct three-dimensional models in the macroscopic or microscopic field. Extracting spatial features for three-dimensional reconstruction from two-dimensional training images (TIs), and characterizing non-stationary features with directional ductility are two key issues in MPS simulation. This study presents a step-wise MPS-based three-dimensional structures reconstruction algorithm with the sequential process and hierarchical strategy based on two-dimensional images. An extension method is proposed to construct three-dimensional TIs. With a sequential simulation process, an initial guess at the coarsest scale is simulated, in which hierarchical strategy is used according to the characteristics of TIs. To obtain a more refined realization, an expectation-maximization like iterative process with global optimization is implemented. A concrete example of chondrite micro-structure simulation, in which one scanning electron microscopy (SEM) image of the Heyetang m eteorite is used as TI, shows that the presented algorithm can simulate complex non-stationary structures.

Multiple-Point Statistics Simulation Models: Pretty Pictures or Decision-Making Tools?

Abundant literature has been produced for the last two decades about multiple-point statistics simulation, or MPS. The idea behind MPS is very simple: reproduce patterns from a 2D, or most often a 3D, training image that displays the type of geological heterogeneity deemed to be relevant to the reservoir or field under study, while honoring local data. Replicating an image is a traditional computer science problem. Thus, it should come as no surprise if a growing number of publications on MPS borrow ideas and techniques directly from computer vision and machine learning to improve the reproduction of training patterns. However, quoting Andre Journel, “Geostatistics is not about generating pretty pictures.” Models have a purpose. For example, in oil and gas applications, reservoir models are used to estimate hydrocarbon volumes and book reserves, run flow simulations to forecast hydrocarbon production and ultimate recovery, and make decisions about field development or optimal well drilling locations. Specific key features such as the extent and connectivity of shale barriers may have a major impact on the reservoir performance forecasts and the field development decisions to be made. Those key features that need to be captured in the model, along with the available subsurface data and constraints of the project, should be the primary drivers in selecting the most appropriate modeling techniques and options to obtain reliable results and make sound decisions. In this paper, the practitioners’ point of view is used to evaluate alternative MPS implementations and highlight remaining gaps.

3D multiple-point statistics simulations of the Roussillon Continental Pliocene aquifer using DeeSse

Abstract. This study introduces a novel workflow to model the heterogeneity of complex aquifers using the multiple-point statistics algorithm DeeSse. We illustrate the approach by modeling the Continental Pliocene layer of the Roussillon aquifer in the region of Perpignan (southern France). When few direct observations are available, statistical inference from field data is difficult if not impossible and traditional geostatistical approaches cannot be applied directly. By contrast, multiple-point statistics simulations can rely on one or several alternative conceptual geological models provided using training images (TIs). But since the spatial arrangement of geological structures is often non-stationary and complex, there is a need for methods that allow to describe and account for the non-stationarity in a simple but efficient manner. The main aim of this paper is therefore to propose a workflow, based on the direct sampling algorithm DeeSse, for these situations. The conceptual model is provided by the geologist as a 2D non-stationary training image in map view displaying the possible organization of the geological structures and their spatial evolution. To control the non-stationarity, a 3D trend map is obtained by solving numerically the diffusivity equation as a proxy to describe the spatial evolution of the sedimentary patterns, from the sources of the sediments to the outlet of the system. A 3D continuous rotation map is estimated from inferred paleo-orientations of the fluvial system. Both trend and orientation maps are derived from geological insights gathered from outcrops and general knowledge of processes occurring in these types of sedimentary environments. Finally, the 3D model is obtained by stacking 2D simulations following the paleo-topography of the aquifer. The vertical facies transition between successive 2D simulations is controlled partly by the borehole data used for conditioning and by a sampling strategy. This strategy accounts for vertical probability of transitions, which are derived from the borehole observations, and works by simulating a set of conditional data points from one layer to the next. This process allows us to bypass the creation of a 3D training image, which may be cumbersome, while honoring the observed vertical continuity.

Quantitative Analysis for the Reconstruction of Porous Media Using Multiple-Point Statistics

The pore structure reconstruction of the porous media is of great importance to the research of mechanisms of fluid flow in porous media. To capture the large-scale patterns in the pore space, the multiple-point statistical technique is generally adopted for porous media reconstruction. Commonly, two different schemes, i.e., the single-grid scheme and the multiple-grid scheme, can be applied for simulation realization. The selection between this two schemes and a proper data template size have thus become a new research issue, and the performance of the characteristic reproduction of the training image using this two schemes must be quantified. In this paper, a series of multiple-point statistics simulation basing on a 2D micro-CT sandstone image are proceeded using both single- and multiple-grid schemes, and different data templates are adapted for porous media reconstruction. Further, to quantify the impact of the computational schemes and setting of the data template to the simulation realizations, a number of measurements considering the pore diameter, porosity, connectivity, and permeability are implemented to fully analyze the results obtained. Results show that by using the single-point statistical method, a large template is necessary to reproduce large-scale structures. The multiple-grid template method may bring great benefits to simulation efficiency over the simple data template method, as well as the recovery of the pore long-range geometric features and seepage characteristics. With the extension of the template for the multiple-grid scheme, the simulation results show lack of variations to some extent.

Recursive convolutional neural networks in a multiple-point statistics framework

This work proposes a new technique for multiple-point statistics simulation based on a recursive convolutional neural network approach coined RCNN. The work focuses on methodology and implementation rather than performance to demonstrate the potential of deep learning techniques in geosciences. Two and three dimensional case studies are carried out. A sensitivity analysis is presented over the main RCNN structural parameters using a well-known training image of channel structures in two dimensions. The optimum parameters found are applied into image reconstruction problems using two other training images. A three dimensional case is shown using a synthetic lithological surface-based model. The quality of realizations is measured by statistical, spatial and accuracy metrics. The RCNN method is compared to standard MPS techniques and an improving framework is proposed by using the RCNNE-type as secondary information. Strengths and weaknesses of the methodology are discussed by reviewing the theoretical and practical aspects.

Multiple-point statistics using multi-resolution images

Multiple-point statistics (MPS) is a simulation technique allowing to generate images that reproduce the spatial features present in a training image (TI). MPS algorithms consist in sequentially filling a simulation grid such that patterns around the simulated values come from the TI. Following this principle, joint simulations of multiple variables can be handled and complex heterogeneous fields can be generated. However, inconsistent patterns are often found in the results and some spatial features can be difficult to reproduce. In this paper, a new MPS algorithm based on a multi-resolution representation of the TI is proposed to enhance the quality of the realizations. The method consists in first building a pyramid of images from the TI by successive convolution using Gaussian-like kernels. Secondly, a MPS simulation is done at the lowest resolution level. Then, the result is expanded to the next level of resolution (one rank higher) and used as a conditioning variable for a joint MPS simulation at that level. This last step is repeated up to the initial resolution, where the final simulation is retrieved. The method is implemented in the DeeSse code based on the direct sampling algorithm. Most of the features provided by the direct sampling (conditioning to hard data, uni- or multi-variate simulation of categorical and continuous variables, scaling and rotation of the training structures) are compatible with the proposed method and the usability is maintained. Finally, various examples show that in most of the situations, combining Gaussian pyramids with MPS allows to get results of better quality and in less time compared to direct MPS simulations.

Multiresolution Approach to Condition Categorical Multiple‐Point Realizations to Dynamic Data With Iterative Ensemble Smoothing

AbstractA new methodology is presented for the conditioning of categorical multiple‐point statistics (MPS) simulations to dynamic data with an iterative ensemble smoother (ES‐MDA). The methodology relies on a novel multiresolution parameterization of the categorical MPS simulation. The ensemble of latent parameters is initially defined on the basis of the coarsest‐resolution simulations of an ensemble of multiresolution MPS simulations. Because this ensemble is non‐multi‐Gaussian, additional steps prior to the computation of the first update are proposed. In particular, the parameters are updated at predefined locations at the coarsest scale and integrated as hard data to generate a new multiresolution MPS simulation. The performance of the methodology was assessed on a synthetic groundwater flow problem inspired from a real situation. The results illustrate that the method converges towards a set of final categorical realizations that are consistent with the initial categorical ensemble. The convergence is reliable in the sense that it is fully controlled by the integration of the ES‐MDA update into the new conditional multiresolution MPS simulations. Thanks to a massively reduced number of parameters compared to the size of the categorical simulation, the identification of the geological structures during the data assimilation is particularly efficient for this example. The comparison between the estimated uncertainty and a reference estimate obtained with a Monte Carlo method shows that the uncertainty is not severely reduced during the assimilation as is often the case. The connectivity is successfully reproduced during the iterative procedure despite the rather large distance between the observation points.

3D stochastic modeling framework for Quaternary sediments using multiple-point statistics: A case study in Minjiang Estuary area, southeast China

Multiple-point statistics (MPS) has shown promise in representing complicated subsurface structures, such as the sedimentary system. The Quaternary sedimentary system possesses prominent anisotropy characteristics. For a practical three-dimensional (3D) application, therefore, the critical challenges come not only from the difficulty of obtaining credible 3D training images, but also from the serious non-stationarity. Moreover, MPS-based simulations are usually performed on a regular Cartesian grid so that they cannot realistically reflect the actual shape and distribution of geological structures. In our work, an integrated MPS-based 3D modeling framework is presented by incorporating the characteristics of Quaternary sediments and the datasets obtained from geological exploration. The framework is embedded into a 3D modeling grid to achieve more precise visualization for the subsurface structures. Following the proposed workflow, we perform the 3D modeling practices for Quaternary sedimentary facies and Quaternary stratigraphic structures by using the datasets from a Quaternary geological survey project at a southeast coastal city in China. The multiple-scale models are constructed by performing the MPS simulations which take into account the constraints of structural conditions. The application verifies the rationality and applicability of the presented workflow where the various structural features are integrated into a unified 3D model to achieve a more realistic characterization.

Correction to: Downscaling Images with Trends Using Multiple-Point Statistics Simulation: An Application to Digital Elevation Models

The original version of this article unfortunately contained mistakes in equations 9, 10, 12 and 13.

Downscaling Images with Trends Using Multiple-Point Statistics Simulation: An Application to Digital Elevation Models

Remote sensing and geophysical imaging techniques are often limited in terms of spatial resolution. This prevents the characterization of physical properties and processes at scales finer than the spatial resolution provided by the imaging sensor. In the last decade, multiple-point statistics simulation has been successfully used for downscaling problems. In this approach, the missing fine-scale structures are imported from a training image which describes the correspondence between coarse and equivalent fine-scale structures. However, in many cases, large variations in the amplitude of the imaged physical attribute, known as trends, pose a challenge for the detection and simulation of these fine-scale features. Here, we develop a novel multiple-point statistics simulation method for downscaling coarse-resolution images with trends. The proposed algorithm relies on a multi-scale sequential simulation framework. Trends in the data are handled by an inbuilt decomposition of the target variable into a deterministic trend component and a stochastic residual component at multiple scales. We also introduce the application of kernel weighting for computing distances between data events and probability aggregation operations for integrating different support data based on a distance-to-probability transformation function. The algorithm is benchmarked against two-point and multiple-point statistics simulation methods, and a deterministic interpolation technique. Results show that the approach is able to cope with non-stationary data sets and scenarios in which the statistics of the training image differ from the conditioning data statistics. Two case studies using digital elevation models of mountain ranges in Switzerland illustrate the method.

Combining Clustering Methods With MPS to Estimate Structural Uncertainty for Hydrological Models

This study presents a novel expansion of the clay-fraction (CF)/resistivity clustering method aiming at developing realizations of subsurface structures based on multiple point statistics (MPS). The CF-resistivity clustering method is used to define a data driven training image (TI) for MPS simulations. By combining this TI with uncertainty estimates obtained from correlation between the resistivity models and the unique categories in the TI, subsurface realizations are generated honoring all geophysical and lithological data including any equivalent resistivity layer sequences. The generated subsurface realizations were calibrated in a steady state groundwater model. Forecasts of well catchment zones were derived based on two wells located in areas with different levels of structural uncertainty. The catchment probability maps derived from the structural realizations were compared with the well catchment forecasted by a deterministic subsurface structure, and we are able to capture this catchment within the estimated uncertainties. We believe that this study is the first to combine MPS methods with a complete data driven workflow going directly from lithological and geophysical data to realizations of the subsurface structures. The main benefits of this is that it is data driven, fast, reproducible, and transparent.

Multiple-point statistical simulation of the ore boundaries for a lateritic bauxite deposit

Resource estimation of mineral deposits requires spatial modelling of orebody boundaries based on a set of exploration borehole data. Given lateritic bauxite deposits, the spacing between the boreholes is often determined based on the grade continuity. As a result, the selected drill spacing might not capture the underlying (true) lateral variability apparent in the orebody boundaries. The purpose of this study is to investigate and address the limitations imposed by such problems in lateritic metal deposits through multiple-point statistics (MPS) framework. Rather than relying on a semivariogram model, we obtain the required structural information from the footwall topographies exposed after previous mining operations. The investigation utilising the MPS was carried out using the Direct Sampling (DS) MPS algorithm. Two historical mine areas along with their mined-out surfaces and ground penetrating radar surveys were incorporated as a bivariate training image to perform the MPS simulations. In addition, geostatistical simulations using the Turning Bands method were also performed to make the comparison against the MPS results. The performances were assessed using several statistical indicators including higher-order spatial cumulants. The results have shown that the DS can satisfactorily simulate the orebody boundaries by using prior information from the previously mined-out areas.

Analogue GPR study of the Permian fanglomerates from Zygmuntówka Quarry near Chęciny, Holy Cross Mountains, southern Poland, for construction of a training image for multiple point simulations

The distribution of Permian alluvial fan lithofacies in a quarry at Zygmuntowka near Checiny, Holy Cross Mts., in southern Poland was investigated using ground penetrating radar (GPR) in order to create a training image for multiple point statistics (MPS) reconstructions of alluvial fan sedimentary facies. Five pseudo-3D GPR datasets were collected, processed and uploaded for interpretation into SKUA-GOCAD 3D geological modelling software. Three radar facies were distinguished based on the 3D geometrical pattern of radar reflections and linked to lithofacies described from the quarry by Zbroja et al. (1998). A statistical summary showed that ~50% of the lithofacies resulted from gravity flows (mostly non-cohesive), while the remaining proportion was deposited by unconfined and confined flash floods. Fluvial sedimentary facies left by waning of catastrophic floods or reworking during fair weather, although not prevalent, could not be distinguished from confined flood deposits based only on GPR data. The GPR datasets together with information from field observations were used to carry out MPS simulations and estimate the most probable 3D model of lithofacies at the quarry scale. This model will in turn serve as a training image for MPS reconstructions of alluvial-fan facies of Rotliegend conglomerates in the multi-scale geological model of the Gorzow Block (western Poland).

Enhanced Multiple-Point Statistical Simulation with Backtracking, Forward Checking and Conflict-Directed Backjumping

During a conventional multiple-point statistics simulation, the algorithm may not find a matched neighborhood in the training image for some unsimulated pixels. These pixels are referred to as the dead-end pixels; the existence of the dead-end pixels means that multiple-point statistics simulation is not a simple sequential simulation. In this paper, the multiple-point statistics simulation is cast as a combinatorial optimization problem, and the efficient backtracking algorithm is developed to solve this optimization problem. The efficient backtracking consists of backtracking, forward checking, and conflict-directed backjumping algorithms that are introduced and discussed in this paper. This algorithm is applied to simulate multiple-point statistics properties of some synthetic training images; the results show that no anomalies occurred in any of the produced realizations as opposed to previously published methods for solving the dead-end pixels. In particular, in simulating a channel system, all the channels generated by this method are continuous, which is of paramount importance in fluid flow simulation applications. The results also show that the presence of hard data does not degrade the quality of the generated realizations. The presented method provides a robust algorithmic framework for performing MPS simulation.

Improving Filtersim Simulation Algorithm of Multiple-Point Geostatistics Using New Training Image Based Adaptive Filters and Supervised Hard Data Conditioning

Summary Simulations based on Multiple-point statistics (MPS) form a class of newly developed techniques that have received great attention in recent years and are mainly used to map spatial complexity and heterogeneity of geological phenomenon. One of the most commonly used pattern based multiple-point geostatistical simulation algorithm is called Filtersim. In the conventional Filtersim algorithm, the detected patterns in training images are transformed into filter score space using a fixed number of filters that are neither dependent on the training images nor on the characteristics of patterns extracted from them. Through using principal component analysis in current study, a set of new filters are designed in such a way to include most structural information specific to the governing training image resulting in the selection of closer patterns in the filter score space. Comparing the results of applying our proposed filters to that of conventional Filtersim algorithm shows a significant improvement in recovering expected shapes and structural continuity in the final simulated realizations. Introduction MPS simulation is capable of reproducing complex geological patterns that cannot be modeled by two-point statistics moments such as variograms. One of the commonly used pattern-based approach, FILTERSIM, introduces filters to summarize high dimensional patterns into a filter score space. Then, it classifies the patterns into a limited number of classes. FILTERSIM uses six predefined filters, which are not specific to the Training Image. The aim of current research is to introduce new improvements to the initial FILTERSIM technique through defining new adaptive filters. Methodology and Approaches In contrast to original FILTERSIM method, our approach uses some filters which are designed specifically for any given Training Image using PCA analysis. In addition, our approach uses a combination of raster path and random partitioning methods in the course of simulation. The raster path approach results in simulations showing a good continuity of the patterns. The random partitioning is a new feature that allows better simulated realizations as compared to the straightforward use of the raster path. Results and Conclusions The newly proposed approach have been tested on some training images and the results have been compared with previous pattern-based algorithms. Through visual inspection and some discrepancy measures it is shown that the simulated realizations obtained by our approach are much closer to the training image than those obtained with the other methods. In particular, the continuity of the patterns is better preserved with the proposed method.