Sparsity Constraint

As a linear transformation, the apex-shifted Radon transform can be represented simply as  

$$\mathbf{d}=\mathbf{Lm} ,$$ (1)

where ${\bf d}$ is the image in the angle domain, ${\bf m}$ is the image in the Radon domain and ${\bf L}$ is the forward apex-shifted Radon transform operator. To find the model ${\bf m}$ that best fits the data in a least-squares sense, we minimize the objective function:  

$$f({\bf m}) = \Vert{\bf Lm-d}\Vert^2 + \epsilon^2b^2 \sum_{i=1}^n \ln\Big(1 +\frac{m_i^2}{b^2}\Big),$$ (2)

where the second term is a Cauchy regularization that enforces sparseness in the model space. Here n is the size of the model space, and ${\epsilon}$ and b are two constants chosen a-priori: ${\epsilon}$ which controls the amount of sparseness in the model space and b related to the minimum value below which everything in the Radon domain should be zeroed Sava and Guitton (2003). The least-squares inverse of ${\bf m}$ is

$${\bf \hat{m}} = \left [ {\bf L'L}+\epsilon^2 {\bf diag}\Big(\frac{1}{1+\frac{m_i^2}{b^2}}\Big) \right ]^{-1}{\bf L'd},$$ (3)

where ${\bf diag}$ defines a diagonal operator. Because the model space can be large, we estimate ${\bf m}$ iteratively. Notice that the objective function in Equation (2) is non-linear because the model appears in the definition of the regularization term. Therefore, we use a limited-memory quasi-Newton method Guitton and Symes (2003) to find the minimum of $f({\bf m})$.