Simulating a depositional (or stratigraphic) sequence conditionally on borehole data is a long-standing problem in hydrogeology and in petroleum geostatistics. This paper presents a new rule-based approach for simulating depositional sequences of surfaces conditionally on lithofacies thickness data. The thickness of each layer is modeled by a transformed latent Gaussian random field allowing for null thickness thanks to a truncation process. Layers are sequentially stacked above each other following the regional stratigraphic sequence. By choosing adequately the variograms of these random fields, the simulated surfaces separating two layers can be continuous and smooth. Borehole information is often incomplete in the sense that it does not provide direct information about the exact layer that some observed thickness belongs to. The latent Gaussian model proposed in this paper offers a natural solution to this problem by means of a Bayesian setting with a Markov chain Monte Carlo (MCMC) algorithm that can explore all possible configurations that are compatible with the data. The model and the associated MCMC algorithm are validated on synthetic data and then applied to a subsoil in the Venetian Plain with a moderately dense network of cored boreholes.