Evaluating $\int_0^\infty \left| \frac{\sin t}{t} \right|^n \, \mathrm{d}t$ for $n = 3, 5, 7, \dots$

I would like to determine the general term of the following sequence defined by an infinite integral: $$ I_n = \int_0^\infty \left| \frac{\sin t}{t} \right|^n \, \mathrm{d}t \, , $$ wherein $n =3, 5, 7, \dots$ is an odd integer.

It can be checked that the integral is convergent for all values of $n$ in the prescribed range. The case of even $n$ is solved in A sine integral $\int_0^{\infty} \left(\frac{\sin x }{x }\right)^n\,\mathrm{d}x$. Also, $I_1 = \infty$. I have tried to use the method of multiple integrations by parts but in vein. I was wondering whether there exists a suitable approach to address this problem more effectively.

Here is a partial answer: If $n \geq 2$ is an integer, then

$$ I_n = \frac{1}{(n-1)!2^{n-1}} \sum_{l=0}^{\lfloor n/2 \rfloor} (-1)^l \binom{n}{l} (n-2l)^{n-1} J_{n-2l}, \tag{1} $$

where $J_p$ is defined by

$$ J_{p} = \begin{cases} \displaystyle \frac{4p}{\pi} \sum_{j=1}^{\infty} \frac{\log(2\pi j)}{4j^2-p^2}, & \text{if $p$ is odd}, \\ \displaystyle \frac{\pi}{2}, & \text{if $p$ is even}. \end{cases} $$

Proof of $\text{(1)}$. The case of even $n$ has already been discussed in other postings, so we focus on odd $n$. We first note that, if $n \geq 1$ is an odd integer, then

\begin{align*} \frac{\mathrm{d}^{n-1}}{\mathrm{d}x^{n-1}} \sin^n x &= \frac{1}{(2i)^n} \sum_{l=0}^{\frac{n-1}{2}} (-1)^l \binom{n}{l} \frac{\mathrm{d}^{n-1}}{\mathrm{d}x^{n-1}} (e^{(n-2l)ix} - e^{-(n-2l)ix}) \\ &= \frac{1}{2^{n-1}} \sum_{l=0}^{\frac{n-1}{2}} (-1)^l \binom{n}{l} (n-2l)^{n-1} \sin((n-2l)x). \end{align*}

So by applying integration by parts $(n-1)$-times, we get

\begin{align*} I_n &= \sum_{k=0}^{\infty} \int_{0}^{\pi} \frac{\sin^n x}{(x+k\pi)^n} \, \mathrm{d}x \\ &= \frac{1}{(n-1)!} \sum_{k=0}^{\infty} \int_{0}^{\pi} \biggl( \frac{1}{x+k\pi} - \frac{1}{(k+1)\pi} \biggr) \biggl( \frac{\mathrm{d}^{n-1}}{\mathrm{d}x^{n-1}} \sin^n x \biggr) \, \mathrm{d}x \\ &= \frac{1}{(n-1)!2^{n-1}} \sum_{l=0}^{\frac{n-1}{2}} (-1)^l \binom{n}{l} (n-2l)^{n-1} J_{n-2l}, \end{align*}

where $J_{p}$ is defined by

\begin{align*} J_{p} &= \sum_{k=0}^{\infty} \int_{0}^{\pi} \biggl( \frac{1}{x+k\pi} - \frac{1}{(k+1)\pi} \biggr) \sin(px) \, \mathrm{d}x. \end{align*}

If $p$ is odd, then the above definition is recast as

\begin{align*} J_{p} &= \sum_{k=0}^{\infty} \biggl( \int_{0}^{\pi} \frac{1}{x+k\pi} \sin(px) \, \mathrm{d}x - \frac{2}{p\pi(k+1)} \biggr) \\ &= \lim_{N \to \infty} \biggl( \int_{0}^{N \pi} \frac{\sin(p(x \text{ mod } \pi))}{x} \, \mathrm{d}x - \frac{2}{p\pi} H_N \biggr), \end{align*}

where $H_N = 1 + \frac{1}{2} + \dots + \frac{1}{N}$ is the $N$-th harmonic number. Still assuming that $p$ is an odd integer, Fourier series computation shows that

\begin{align*} \sin(p(x \text{ mod } \pi)) &= \frac{2}{p\pi} - \frac{4p}{\pi} \sum_{n=1}^{\infty} \frac{\cos(2\pi n x)}{4n^2-p^2} \\ &= \frac{4p}{\pi} \sum_{j=1}^{\infty} \frac{1 - \cos(2\pi j x)}{4j^2-p^2}, \end{align*}

and so,

\begin{align*} \int_{0}^{N \pi} \frac{\sin(p(x \mathrm{ mod } \pi))}{x} \, \mathrm{d}x &= \frac{4p}{\pi} \sum_{j=1}^{\infty} \frac{1}{4j^2-p^2} \int_{0}^{N \pi} \frac{1 - \cos(2 j x)}{x} \, \mathrm{d}x \\ &= \frac{4p}{\pi} \sum_{j=1}^{\infty} \frac{1}{4j^2-p^2} (\gamma + \log(2\pi j N) - \operatorname{Ci}(2\pi j N) ). \end{align*}

Plugging this back and using the identity $\frac{4p}{\pi} \sum_{j=1}^{\infty} \frac{1}{4j^2-p^2} = \frac{2}{p\pi}$, which itself follows from the Fourier series of $\sin(p(x \text{ mod } \pi))$, we finally obtain

\begin{align*} J_{p} &= \frac{4p}{\pi} \lim_{N \to \infty} \sum_{j=1}^{\infty} \frac{1}{4j^2-p^2} (\gamma + \log(2\pi j N) - \operatorname{Ci}(2\pi j N) - H_N ) \biggr) \\ &= \frac{4p}{\pi} \sum_{j=1}^{\infty} \frac{\log(2\pi j)}{4j^2-p^2} \end{align*}

as desired. $\square$

Addendum. Here is a Mathematica code for numerical verification of $\text{(1)}$:

n = 5; (* Choose your favorite odd integer >= 3*)
NIntegrate[Evaluate[Sum[1/(x + k Pi)^n, {k, 0, Infinity}] Sin[x]^n], {x, 0, Pi}, WorkingPrecision -> 20]
TermJ[p_] := (4 p)/Pi NSum[Log[2 Pi j]/(4 j^2 - p^2), {j, 1, Infinity}, WorkingPrecision -> 20];
1/((n - 1)! 2^(n - 1)) Sum[Binomial[n, l] (-1)^l (n - 2 l)^(n - 1) TermJ[n - 2 l], {l, 0, (n - 1)/2}]
Clear[n, TermJ];