Marcus Carlsson

Viktor V. Nikitin, Fredrik Andersson, Marcus Carlsson et al.

The hyperbolic Radon transform is a commonly used tool in seismic processing, for instance in seismic velocity analysis, data interpolation and for multiple removal. A direct implementation by summation of traces with different moveouts is computationally expensive for large data sets. In this paper we present a new method for fast computation of the hyperbolic Radon transforms. It is based on using a log-polar sampling with which the main computational parts reduce to computing convolutions. This allows for fast implementations by means of FFT. In addition to the FFT operations, interpolation procedures are required for switching between coordinates in the time-offset; Radon; and log-polar domains. Graphical Processor Units (GPUs) are suitable to use as a computational platform for this purpose, due to the hardware supported interpolation routines as well as optimized routines for FFT. Performance tests show large speed-ups of the proposed algorithm. Hence, it is suitable to use in iterative methods, and we provide examples for data interpolation and multiple removal using this approach.

Fredrik Andersson, Marcus Carlsson

We present a new method for complex frequency estimation in several variables, extending the classical (1d) ESPRIT-algorithm. We also consider how to work with data sampled on non-standard domains (i.e going beyond multi-rectangles).

Fredrik Andersson, Marcus Carlsson

We develop fixed-point algorithms for the approximation of structured matrices with rank penalties. In particular we use these fixed-point algorithms for making approximations by sums of exponentials, or frequency estimation. For the basic formulation of the fixed-point algorithm we show that it converges to the minimum of the convex envelope of the original objective function along with its structured matrix constraint. It often happens that this solution agrees with the solution to the original minimization problem, and we provide a simple criterium for when this is true. We also provide more general fixed-point algorithms that can be used to treat the problems of making weighted approximations by sums of exponentials given equally or unequally spaced sampling. We apply the method to the case of missing data, although optimal convergence of the fixed-point algorithm is not guaranteed in this case. However, it turns out that the method often gives perfect reconstruction (up to machine precision) in such cases. We also discuss multidimensional extensions, and illustrate how the proposed algorithms can be used to recover sums of exponentials in several variables, but when samples are available only along a curve.

Fredrik Andersson, Marcus Carlsson

We consider sequences $(B_k)_{k=0}^\infty$ of points obtained by projecting back and forth between two manifolds $\M_1$ and $\M_2$, and give conditions guaranteeing that the sequence converge to a limit $B_\infty\in\M_1\cap\M_2$. Our motivation is the study of algorithms based on finding the limit of such sequences, which have proven useful in a number of areas. The intersection is typically a set with desirable properties, but for which there is no efficient method of finding the closest point $B_{opt}$ in $\M_1\cap\M_2$. We prove not only that the sequence of alternating projections converges, but that the limit point is fairly close to $B_{opt}$, in a manner relative to the distance $\|B_0-B_{opt}\|$, thereby significantly improving earlier results in the field. A concrete example with applications to frequency estimation of signals is also presented.

Fredrik Andersson, Marcus Carlsson, Per-Anders Ivert

The problem of approximating a sampled function using sums of a fixed number of complex exponentials is considered. We use alternating projections between fixed rank matrices and Hankel matrices to obtain such an approximation. Convergence, convergence rates and error estimates for this technique are proven, and fast algorithms are developed. We compare the numerical results obtain with the MUSIC and ESPRIT methods.

Fredrik Andersson, Marcus Carlsson, Viktor V. Nikitin

The Radon transform and its adjoint, the back-projection operator, can both be expressed as convolutions in log-polar coordinates. Hence, fast algorithms for the application of the operators can be constructed by using FFT, if data is resampled at log-polar coordinates. Radon data is typically measured on an equally spaced grid in polar coordinates, and reconstructions are represented (as images) in Cartesian coordinates. Therefore, in addition to FFT, several steps of interpolation have to be conducted in order to apply the Radon transform and the back-projection operator by means of convolutions. Both the interpolation and the FFT operations can be efficiently implemented on Graphical Processor Units (GPUs). For the interpolation, it is possible to make use of the fact that linear interpolation is hard-wired on GPUs, meaning that it has the same computational cost as direct memory access. Cubic order interpolation schemes can be constructed by combining linear interpolation steps which provides important computation speedup. We provide details about how the Radon transform and the back-projection can be implemented efficiently as convolution operators on GPUs. For large data sizes, speedups of about 10 times are obtained in relation to the computational times of other software packages based on GPU implementations of the Radon transform and the back-projection operator. Moreover, speedups of more than a 1000 times are obtained against the CPU-implementations provided in the MATLAB image processing toolbox.