  We propose a sparse regularization model for inversion of incomplete Fourier
transforms and apply it to seismic wavefield modeling. The objective function
of the proposed model employs the Moreau envelope of the $\ell_0$ norm under a
tight framelet system as a regularization to promote sparsity. This model leads
to a non-smooth, non-convex optimization problem for which traditional
iteration schemes are inefficient or even divergent. By exploiting special
structures of the $\ell_0$ norm, we identify a local minimizer of the proposed
non-convex optimization problem with a global minimizer of a convex
optimization problem, which provides us insights for the development of
efficient and convergence guaranteed algorithms to solve it. We characterize
the solution of the regularization model in terms of a fixed-point of a map
defined by the proximity operator of the $\ell_0$ norm and develop a
fixed-point iteration algorithm to solve it. By connecting the map with an
$\alpha$-averaged nonexpansive operator, we prove that the sequence generated
by the proposed fixed-point proximity algorithm converges to a local minimizer
of the proposed model. Our numerical examples confirm that the proposed model
outperforms significantly the existing model based on the $\ell_1$-norm. The
seismic wavefield modeling in the frequency domain requires solving a series of
the Helmholtz equation with large wave numbers, which is a computationally
intensive task. Applying the proposed sparse regularization model to the
seismic wavefield modeling requires data of only a few low frequencies,
avoiding solving the Helmholtz equation with large wave numbers. Numerical
results show that the proposed method performs better than the existing method
based on the $\ell_1$ norm in terms of the SNR values and visual quality of the
restored synthetic seismograms.
