This paper presents a fast algorithm to solve a spectral estimation problem for two-dimensional random fields. The latter is formulated as a convex optimization problem with the Itakura-Saito pseudodistance as the objective function subject to the constraints of moment equations. We exploit the structure of the Hessian of the dual objective function in order to make possible a fast Newton solver. Then we incorporate the Newton solver to a predictor-corrector numerical continuation method which is able to produce a parametrized family of solutions to the moment equations. We have performed two sets of numerical simulations to test our algorithm and spectral estimator. The simulations on the frequency estimation problem shows that our spectral estimator outperforms the classical windowed periodograms in the case of two hidden frequencies and has a higher resolution. The other set of simulations on system identification indicates that the numerical continuation method is more robust than Newton's method alone in ill-conditioned instances.