In this work we present a novel level set technique for shape reconstruction in history matching for reservoirs with two or more kinds of rocks (the so-called lithofacies) using stochastic initializations. In the paper we discuss the use of sequential Gaussian simulation for the creation of geostatistical initial guesses which will then be applied to an earlier introduced level set based shape reconstruction algorithm for the history matching problem in reservoir characterization. The shapes or regions with sharp interfaces to each other are represented by a level set function which needs to be determined from the production data gathered at few locations of the reservoir. We present numerical results in 2D which demonstrate that our method is able to provide reliable estimates of this structure from these relatively few production data even though the topology of the unknown regions is a priori unknown. © 2008 Elsevier Ltd. All rights reserved.