How to derive the Hankel transform of the following function $f(r)$?
$$ f(r) = \frac{1}{\sqrt{r^2+a^2}}\,,\quad a\in\mathbb{R}$$
Since Hankel transform is an involution, it is equivalent to show that the Hankel transform of $g(k)$ (I know the answer) is the $f(x)$ where
$$g(k) = \frac{e^{|k|a}}{|k|}\,$$
(up to some constant). Hankel transform $\mathcal{H}$ of a function $f$ is defined as:
$${ \mathcal{H}[f](k)=\int _{0}^{\infty }f(r)J_{0}(kr)\,r\,\mathrm {d} r}\,,$$
where $J_0$ is a zeroth-order Bessel function.
I know the result because the transform is given in several Hankel transform tables, like here. I also saw some more explanation in the following answer (which gives the transform result in a more general case). However, I am interested in the derivation of this result, could you refer me to a paper where the proof is given or hint how one might approach this integration?
It strikes me that this quite complex integral with the Bessel function form reduces to simple functions like exponential and square root.
The $g$ integral has a nice method. We have $$ \mathcal{H}[g](r) = \int_0^\infty \frac{e^{-ak}}{k}J_0(kr)kdk = \frac{1}{r}\int_0^\infty e^{-(a/r)u}J_0(u)du. $$ This is the Laplace transform of $J_0$. So let's solve the Laplace transformed-Bessel equation. This is particularly convenient for the order $0$ case, as we can cancel out an $x$ from each term. This gives $$ \mathcal{L}\left[x\frac{d^2J_0}{dx^2}+\frac{dJ_0}{dx} +xJ_0\right](s) = -\frac{d}{ds}(s^2\tilde{J}_0) + s\tilde{J}_0 -\frac{d\tilde{J}_0}{ds} = -(s^2+1)\frac{d\tilde{J}_0}{ds}-s\tilde{J}_0 = 0 $$ This separable first order ODE is easily solved to give $\tilde{J}_0 = (1 + s^2)^{-1/2}$. So we have $$ \mathcal{H}[g](r) = \frac{1}{r}\int_0^\infty e^{-(a/r)u}J_0(u)du = \frac{1}{r}\tilde{J}_0\left(\frac{a}{r}\right) = \frac{1}{r\sqrt{1+(a/r)^2}} = \frac{1}{\sqrt{r^2 + a^2}} $$