Although the interpolation of seismic data beyond aliasing is still under development, the theories and techniques for interpolating unaliased data are well understood and routinely used in industry. In the 1950s, Yen (1956 and 1957) developed an algorithm for interpolating or resampling irregularly sampled, one-dimensional data. In his algorithm, the irregularly sampled data are modeled as the interpolation of unknown uniformly sampled data:
Hale (1980) tested Yen's algorithm and pointed out that, although the method is exact for noise-free data, it could give erroneous results for noisy data. With Yen's method, a single noisy sample, even one far away from the correct sampling positions, may have an undesirable influence on the values of the interpolated samples. Moreover, Yen's algorithm is computationally expensive. Resampling N points requires an order of N3 operations. To solve these two problems, Hale replaced the sinc function in equation (1) with a truncated, tapered sinc function. Although this procedure sacrifices accuracy at high frequencies, it performs better and at much less computational cost than Yen's method for data contaminated with truncations, aliasing data or noisy data.
This paper describes a new, one-dimensional interpolation method. In this method, a continuous signal is modeled with its Fourier series expansion. The coefficients of this expansion are discrete spectra of the signal, and can be estimated by fitting the modeled data to the irregularly sampled data using the standard least squares technique. It turns out that the normal matrix in this least squares problem is a Hermitian Toeplitz matrix, and its inverse requires only an order of N2 operations. Therefore, the algorithm is quite efficient. Assuming that the irregularly sampled signal is not aliased and that it has a finite duration, the continuous signal can be recovered from the irregular samples exactly. This algorithm can also be used to compute the discrete Fourier transform of the irregularly sampled data.
In the following sections, I first formulate the problem and observe the Toeplitz structure of the matrix to be inverted. Then, I show two examples of one-dimensional interpolation and one example of computing the discrete Fourier transform. Finally, I offer conclusions and discuss some potential problems of the method.