Accurate subsurface modeling and characterization require the prediction of facies and rock properties within the reservoir model. This is commonly achieved by inverting geophysical data, such as seismic reflection data, using a two-step approach either in the discrete or the continuous domain. We have adopted an iterative simultaneous method, namely, stochastic perturbation optimization, to invert seismic reflection data jointly for facies and rock properties. Facies first are simulated according to a Markov chain model, and then rock properties are generated with stochastic sequential simulation and cosimulation conditioned to each facies. Elastic and seismic data are computed by applying a rock-physics model to the realizations of petrophysical properties and a seismic convolutional model. The similarity between observed and synthetic seismic data is used to update the solution by perturbing facies and rock properties until convergence. Coupling the discrete and continuous domains ensures a consistent perturbation of the reservoir models throughout the iterations. We have evaluated the method in a 1D synthetic example for the estimation of facies and porosity from zero-offset seismic data assuming a linear rock-physics model to demonstrate the validity of the method. Then, we apply the method to a real 3D data set from the North Sea for the joint estimation of facies and petrophysical properties from prestack seismic data. The results show spatially consistent rock and fluid inverted models in which the predicted facies reproduce the vertical ordering as observed in the well data.