Convergence rate of full wavefield inversion using spectral shaping

The invention relates to a method for speeding up iterative inversion of seismic data (106) to obtain a subsurface model (102) by using local cost function optimization. The frequency spectrum of the updated model at each iteration is controlled to match a known or estimated frequency spectrum for the subsurface region, preferably the average amplitude spectrum of the subsurface P-impedance. The controlling is done either by applying a spectral-shaping filter to the source wavelet (303) and to the data (302) or by applying the filter, which may vary with time, to the gradient of the cost function (403). The source wavelet's amplitude spectrum (before filtering) should satisfy D(f)=fI(f)W(f), where f is frequency, D(f) is the average amplitude spectrum of the seismic data, and I(f) is the average amplitude spectrum for P-impedance in the subsurface region (306,402) or an approximation thereof.