# Using the Mellin transform to invert ill-posed problems of harmonic transforms

I will try to explain the problem using just words. Let's see how far I get!

I will use the harmonic transform of the Radon transform for 2 or more dimensions as an example, but there is a large class of problems out there that this might apply to.

The use of the Mellin transform to solve Mellin convolutions of integral equations often results in an inversion that is unstable, in the sense that a numerical discretization of the solution fails to converge or is otherwise unmanageable. As a result, the formal, analytical result is often considered flawed. When this difficulty is encountered, this avenue of attack is usually abandoned. A primary example of this difficulty is the harmonic inversion of the Radon transform in $\mathbb{R}^{n}$.

This inversion is based on the the ultraspherical, spherical, and circular harmonic transforms of the Radon transform that have formed the core of theoretical study in tomography. From these transforms, researchers have been able to determine important information that is utilized in practical and commercial applications of tomography. Cormack's early work (Alan M. Cormack. Journal of Applied Physics, 34(9):2722–2727, 1963.), (ibid. J. Appl. Phys., 35(10):2908–2913, 1964) was critical in demonstrating the feasibility of tomographic reconstruction. Information about sampling rates and the null space of the Radon transform has been obtained from the harmonic form {F.Natterer. Mathematics of Computerized Tomography. John Wiley and Sons, New York, 1986.}, {Natterer and Wubbeling, 2001}, {A. K. Louis. SIAM J. Math. Analy., 15(3):321–633, 1984.}. Practical reconstruction algorithms, though one exists (Hawkins & Barrett, 1986), have been too late arriving on the scene to have any impact on commercial tomographic reconstruction.

Part of the reason for this lack of use was that the inversion of the harmonic solution, in the form of a Mellin convolution integral equation, first obtained for $\mathbb{R}^{n}$ by Louis {A. K. Louis, 1984}, was extremely unstable as an exercise in numerical integration. Numerical integration is often utilized to solve integral equations, especially if the inversion, as an integral transform, is amenable to numerical techniques. This was definitely not the case, and the inverse harmonic solution to the Radon transform was thought, incorrectly, to be as unstable as limited angle tomography. The solutions proposed by Louis (1984) and {S.R. Deans. The Radon Transform and Some of Its Applications. John Wiley and Sons, New York, 1983} were based on this unstable inverse. It was thought for a long time that any formal or analytical solution based upon harmonic analysis must also be unstable. A method of stabilizing the transform was found by Hansen {E. W. Hansen. J. Opt. Soc. Am., 71:304–8, 1981}, but this method exists only for the two dimensional case and has not been generalized to higher dimensions.

Chapman and Carey {C. H. Chapman & P. W. Carey. Inverse Problems, 2:23–49, 1986} derived a stable circular harmonic inverse for the two dimensional case, equivalent to Hansen's, and applied it to an inversion algorithm that produced image quality equivalent to filtered backprojection. The algorithm was a numerical method based on the stable inverse. P. E. Mijnarends {P. E. Mijnarends. Physical Review, 160(3):512–519, 1967} found the stable harmonic inverse for three dimensions. Hawkins, 1983 {W. G. Hawkins. Ph. D. dissertation, U. Arizona, Tucscon, Az, 1983} independently found the same stable harmonic inverse.

Despite the good news that these stable harmonic inverses brought, harmonic analysis is still viewed skeptically in some quarters. It appeared that for the Radon transform, at least, harmonic analysis was at a dead end. Why this issue was not dealt with more urgently is not known.

In 2007, Hawkins {W.G. Hawkins. International Journal of Image and Graphics, 7(1):17–33, 2007} proved that the Tchebyshev/Zernike polynomial solution to the unstable inverse is also a solution to the stable inverse {W. Hawkins. International Journal of Image and Graphics, 7(1):17–33, 2007.}. Therefore, for this polynomial pair, the stable and unstable circular harmonic inverses are equivalent, as indeed they must be; the inverse of a linear integral transform, if it exists, is unique. If they differ, it can be only by a factor that is equivalent to unity. In that article, he again demonstrated that the method of reconstruction with orthogonal polynomials is essentially equivalent to filtered backprojection in both execution speed (numerical complexity) and image quality.

However, there exists this huge difference in numerical stability due to this factor of unity (In Mellin transform space).

Is anyone out there familiar with this phenomenon? Can anyone offer any suggestions? I have found examples of stability/instability for other Mellin transform solutions to ill-posed problems. I also have been able to deduce a general theorem about stability.

-