Intensity inhomogeneity (IIH) in magnetic resonance imaging (MRI) adversely affects image processing and image analysis. In this paper, we propose a novel surface fitting approach for intensity nonuniformity correction in MR images. In our method, pattern matching used in non-local means (NLM) algorithm are used to find voxels with similar structures, and a smooth polynomial surface is fitted by a modified least squares method. This method has the advantage that it can choose appropriate points to easily fit the data. The proposed method has been tested on simulated MR images and real MR brain images, and compared with other popular methods. Both qualitative and quantitative evaluations proved the effectiveness of our proposed method.