In this paper, a numerical simulation method for generating rough surfaces with desired autocorrelation function (ACF) and statistical parameters, including root mean square (rms), skewness (Ssk), and kurtosis (Ku), is developed by combining the polar method, Johnson translation system, and random switching system. This method can be used to generate Gaussian, non-Gaussian, isotropic, and nonisotropic rough surfaces. The simulation results show the excellent performance of present method for producing surface with various desired statistical parameters and ACF. The advantage of present method is that the deviation of statistical parameters and ACF from the desired ones can be as small as required since it is controlled by iterative algorithms.