Consider beautiful integral:-
$I(a, b, c) = \int_0^{\infty} \mathrm{d}x~\sin(a x) I_0(b x) K_1(c x)$
where:
$I_0$ modified Bessel function of first kind (order zero)
$K_1$ modified Bessel function of second kind (order one)
real parameters $a$ and $c > b > 0$
How to evaluate integral?
This is only a partial answer, because I'm not sure how to evaluate the last integral that I arrive at, but I've done some work transforming it into a contour integral that only contains elementary functions. It's also far too long for a comment, so here goes.
We plug in the integral expressions for $I_0(x)$ and $K_1(x)$ and obtain:
\begin{align} I(a,b,c) &= \frac{1}{\pi}\int_0^\infty dx\sin(ax)x\int_1^\infty dy\sqrt{y^2-1}e^{-cxy}\int_0^{\pi}dz e^{bx\cos z}\\ &=\frac{1}{\pi}\Im\left[\int_1^\infty dy\sqrt{y^2-1}\int_0^\pi dz\int_0^\infty dx xe^{(ia + b\cos z - cy)x}\right] \tag1\\ &=\frac{1}{\pi}\Im\left[\int_1^\infty dy\sqrt{y^2-1}\int_0^\pi\frac{dz}{(b\cos z -cy+ia)^2}\right] \tag2\\ \end{align}
The convergence of the integral in $(1)$ is given, because we know that $cy>c>b\geq b\cos z$ for any values of $y$ and $z$ that are going to be considered in the respective integrals.
Let $\gamma := cy - ia$, then $(2)$ becomes
\begin{align} I(b,\gamma)&=\frac{1}{\pi}\Im\left[\int_1^\infty dy\sqrt{y^2-1}\int_0^\pi\frac{dz}{(b\cos z - \gamma)^2}\right] \tag3 \end{align}
Let $u:=\tan\left(\frac{z}{2}\right)$, then $dz = \frac{2}{1+u^2}du$ and $\cos z = \frac{1-u^2}{1+u^2}$. Then, we split the integral into an integral from $0$ to $\frac{\pi}{2}$ and an integral from $\frac{\pi}{2}$ to $\pi$. The bounds on the first integral become $0$ to $\infty$ and the bounds on the second integral become $-\infty$ to $0$. Putting this all together, we have
\begin{align} I(b,\gamma)&=\frac{1}{\pi b^2}\Im\left[\int_1^\infty dy\sqrt{y^2-1}\int_{-\infty}^\infty du\frac{2}{1+u^2}\frac{1}{\left(\frac{1-u^2}{1+u^2}-\frac{\gamma}{b}\right)^2}\right]\\ &=\frac{2}{\pi b}\Im\left[\int_1^\infty dy\sqrt{y^2-1}\int_{-\infty}^\infty \frac{du}{\left[b(1-u^2)-\gamma(1+u^2)\right]^2}\right]\\ &=\frac{2}{\pi b}\Im\left[\int_1^\infty dy\sqrt{y^2-1}\int_{-\infty}^\infty \frac{du}{\left[-u^2(\gamma+b) -(\gamma-b)\right]^2}\right]\\ &=-\frac{2}{\pi b}\Im\left[\int_1^\infty dy\frac{\sqrt{y^2-1}}{(\gamma+b)^2}\int_{-\infty}^\infty \frac{du}{\left[u^2 +\frac{\gamma-b}{\gamma+b}\right]^2}\right] \tag4 \end{align}
Using the well-known result that
$$\int_{-\infty}^\infty\frac{dx}{(x²+a²)^2} = \frac{\pi}{2a^3}$$
for real values of $a$, $(4)$ becomes
\begin{align} I(b,\gamma) &= -\frac{1}{b}\Im\left[\int_1^\infty dy\frac{\sqrt{y^2-1}}{(\gamma+b)^2}\left(\frac{\gamma+b}{\gamma-b}\right)^\frac{3}{2}\right]\\ &=-\frac{1}{b}\Im\left[\int_1^\infty dy \sqrt{\frac{y^2-1}{(\gamma+b)(\gamma-b)^3}}\right] \tag5 \end{align}
Now recall that $\gamma = cy-ia$, thus $y=\frac{\gamma + ia}{c}$, $dy = \frac{d\gamma}{c}$, $\gamma(1)=c-ia$ and $\gamma(\infty) = \infty - ia$. This turns the integral in $(5)$ into a contour integral in the complex plane:
\begin{align} I(a,b,c)&=-\frac{1}{bc^2}\Im\left[\int_{c-ia}^{\infty-ia}d\gamma\sqrt{\frac{(\gamma+ia)^2 - c^2}{(\gamma+b)(\gamma-b)^3}}\right] \tag6 \end{align}
The integrand clearly has poles at $\gamma=\pm b$, so I'm thinking that there is some way to evaluate this using the the Residue Theorem, but I'm not sure how to treat the square root. Perhaps somebody else can take it from here, though.
