Solution using Fourier transform
The integral can be seen as a Fourier transform:
$$
I := \int_0^\infty \frac{\sin^2(x)}{x^2(x^2+1)} \, dx
= \frac12 \int_{-\infty}^{\infty} \frac{\sin(x)}{x} \frac{\sin(x)}{x} \frac{1}{x^2+1} e^{-i0x} \, dx
\\
= \frac12 \mathcal{F}\{ \frac{\sin(x)}{x} \frac{\sin(x)}{x} \frac{1}{x^2+1} \}(0)
.
$$
By the duality of product and convolution under Fourier transform we can rewrite this as
$$
I =
= \frac12 \frac{1}{(2\pi)^2} \left(\mathcal{F}\{ \frac{\sin(x)}{x} \} * \mathcal{F}\{ \frac{\sin(x)}{x} \} * \mathcal{F}\{ \frac{1}{x^2+1} \} \right)(0)
.
$$
Now,
$$
\mathcal{F}\{ \frac{\sin(x)}{x} \} = \pi\mathbf{1}_{[-1,1]}(\xi),
$$
where $\mathbf{1}_A$ is the indicator function on the set $A,$
and
$$
\mathcal{F}\{ \frac{1}{x^2+1} \} = \pi e^{-|\xi|}
.
$$
Thus,
$$
I = \frac{1}{8\pi^2} \left( \pi\mathbf{1}_{[-1,1]} * \pi\mathbf{1}_{[-1,1]} * \pi e^{-|\bullet|} \right)(0)
\\
= \frac{\pi}{8} \left( \mathbf{1}_{[-1,1]} * \mathbf{1}_{[-1,1]} * e^{-|\bullet|} \right)(0)
.
$$
Here,
$$
(\mathbf{1}_{[-1,1]} * \mathbf{1}_{[-1,1]})(x) = (2-|x|) \, \mathbf{1}_{[-2,2]}(x)
$$
so
$$
\left( \mathbf{1}_{[-1,1]} * \mathbf{1}_{[-1,1]} * e^{-|\bullet|} \right)(\xi)
= \int_{-2}^{2} (2-|\eta|) e^{-|\xi-\eta|} \, d\eta
$$
However, we only need to calculate this at $\xi=0$:
$$
I = \left( \mathbf{1}_{[-1,1]} * \mathbf{1}_{[-1,1]} * e^{-|\bullet|} \right)(0)
= \int_{-2}^{2} (2-|\eta|) e^{-|\eta|} \, d\eta
= 2 \int_{0}^{2} (2-\eta) e^{-\eta} \, d\eta = 2 (1 + e^{-2}).
$$
Thus, finally,
$$
I = \frac{\pi}{8}\cdot 2 (1 + e^{-2}) = \frac{\pi}{4}(1 + e^{-2}).
$$