Abstract: Gravity data is always saved and used in grid format, which is also demanded for frequency domain computation as fast fourier transformation. Therefore, gridding method becomes the first and very important step of gravity data processing. Some gridding methods, such as Shepard and Kriging, were introduced into this step and improved the gridding precision. But when it comes to the case of sparse data, the gridding precision were sharply reduced. Moreover, the regular grid densification method was rarely studied. To meet the demands of regular grid densification, the fractal interpolating curve surface theory was introduced and applied in this paper. The computing formulae was deduced in detail. An experiment was carried out and result showed that, compared to the RBF method, the FIC method was better on the precision but worse on the error extremum. It means the FIC method showed bad behavior in data smoothing. Based on the Bouguer anomaly, a “remove-restore” method was applied to solve this problem. Result showed that not only was the precision improved by BAFIC method, but also was the error extremum reduced to a great extent, which made the result more reliable and applicable. The vertical scale parameter which plays a great role in BAFIC method was chosen after times of experiments in this paper, which could be optimized by better and precise methods.