sequences and series - Infinite summation formula for modified Bessel functions of first kind

I was trying to find a closed form for the integral
$$4\int_0^{\pi/2} t \, I_0(2\kappa\cos{t}) dt \; ,$$
where

$$I_{\alpha}(z) := i^{-\alpha}J_{\alpha}(iz) = \sum_{m=0}^{\infty}\frac{\left(\frac{z}{2}\right)^{2m+\alpha}}{m! \Gamma(m+1+\alpha)} = \frac{1}{2\pi} \int_{-\pi}^{\pi} e^{i\alpha \tau + z \sin{\tau}} d\tau$$
are the modified Bessel functions. This integral popped up when I was trying to find the average difference of two points on a circle, where these points are assumed to be drawn independently from a von Mises distribution. It was noted by Robert Israel that this integral can be reduced to

$$\int_0^\pi t I_0(2\kappa \cos(t/2)) \; dt = \frac{\pi^2}{2} I_0(\kappa)^2 - 4 \sum_{r=0}^\infty \frac{I_{2r+1}(\kappa)^2}{(2r+1)^2} \; .$$
So I was wondering, if we can further simplify this expression, or to stated more clearly:

Is there a closed formula for the following sum of modified Bessel functions of the first kind?
$$\sum_{r=0}^\infty \frac{I_{2r+1}(\kappa)^2}{(2r+1)^2}$$

---

A lot of remarkable identities in terms of infinite sums of Bessel functions are known. E.g. Abramowitz and Stegun list in §9.6.33ff. a few of them, like:

$$\begin{align} 1 &= I_0(z) + 2\sum_{r=1}^{\infty} (-1)^{r}I_{2r}(z) \\ e^z &= I_0(z) + 2\sum_{r=1}^{\infty} I_{r}(z) \\ \cosh{z} &= I_0(z) + 2\sum_{r=1}^{\infty} I_{2r}(z) \\ \end{align}$$
WolframResearch lists another bunch of infinite series identities.
Also, Neumann's addition theorem seems to work wonders <1> <2> <3>.

---

Regarding the integral itself, Gradshteyn and Ryzhik mention in 6.519.1 that

$$\int_0^{\pi/2} J_{2r}(2\kappa\cos{t}) = \frac{\pi}{2} J_r^2(\kappa) \; ,$$
where $J_r(x) = i^rI_r(-ix)$. So there might be a chance to expect something along this line.

---

Going back to the original problem "What is the expected value of a distribution $\Delta$ with the following density function"

$$f_{\Delta}(t) := \frac{I_0 \left( 2\kappa \cos{\frac{t}{2}} \right)}{\pi I^2_0(\kappa)} \; ,$$
a straightforward integration leads to the integral mentioned above. Using some probability theory voodoo we can make use of the fact that

$$\mathbb{E}[\Delta] = -i \varphi'_{\Delta}(0) = -i \left[\frac{d}{d\omega} \mathcal{F}(f_{\Delta})(\omega) \right] \Bigg|_{\omega=0} = -i \left[\frac{d}{d\omega} \int_{-\infty}^{\infty} e^{it\omega}f_{\Delta}(t) dt \right] \Bigg|_{\omega=0}$$

where $\varphi_{\Delta}$ is the characteristic function of $f_{\Delta}$ and $\mathcal{F}$ the (properly scaled) Fourier transform. Now with $\varphi(-\omega) = \overline{\varphi(\omega)}$, we could further rewrite

$$\mathbb{E}[\Delta] = -i\varphi'_{\Delta}(0) = \lim_{\omega \rightarrow 0} \frac{\varphi_{\Delta}(\omega) - \varphi_{\Delta}(-\omega)}{2i\omega} = \lim_{\omega \rightarrow 0} \frac{\mathcal{Im}\left(\varphi_{\Delta}(\omega)\right)}{\omega} \,$$

to (by plugging in the integral representation of $I_0$) obtain

$$\mathbb{E}[\Delta] = \frac{\pi}{2} - \frac{4}{\pi I_0^2(\kappa)} \sum_{r=0}^\infty \left( \frac{I_{2r+1}(\kappa)^2}{2r+1} \right)^2 = \frac{1}{\pi^2 I_0^2{\kappa}} \cdot \lim_{\omega \rightarrow 0} \int_0^{\pi/2} \int_{-\pi}^{\pi} \frac{\sin(t\omega)}{\omega} e^{2\kappa\cos{t}\sin{\tau}} d\tau \, dt \; ,$$

but this will essentially lead to the same integral we started with. The promising part about this approach is that the Fourier transform pops up, which might leave some room for the harmonic analysis people among you to do your magic.