Purpose: The purpose of this study is to correct intensity inhomogeneity in MRI. A novel retrospective approach, based on a parametric surface-fitting model, to correct for the smoothly varying inhomogeneity gradient field is described. This method has been evaluated on multiple MR images and has shown promising results. Method and Materials: The major improvement in this study is the calculation of the gradient partial derivatives. Consider the following model of acquired MRI: v(x, y) = g(x, y)u(x, y) + n(x, y), (1) where v(x,y), g(x,y), u(x,y), and n(x,y) are the measured signal, gradient field, uncorrupted signal, and noise, respectively. h+1 samples along the y direction are used to calculate the sum of differences xv in x at a fixed x: (2) Similarly, the sum of h+1 adjacent pixel sums Σxv is: (3) The estimated partial derivative in respect to x is obtained by dividing Eq.(2) by Eq.(3): (4) Partial derivatives are calculated along both x and y directions. The estimated gradient is obtained by integrating the partial derivatives with energy minimization and scale optimization. Finally the estimated gradient is fit to a second order polynomial. Results: The method was tested on checkerboard images with known gradient, as well as MR images of rat brain and spinal cord. For the checkerboard images, restored images were comparable to the uncorrupted images. For MR images, the restored images showed qualitative visual improvement. Conclusions: A novel surface-fitting method was developed to calculate the smooth varying gradient field. The proposed method did not require any a priori knowledge about the gradient nor the anatomy. The corrected results were visually more similar to their uncorrupted forms. Combined, these corrections can improve segmentation and registration quality.