How to evaluate $\displaystyle \int_{0}^{\infty }\frac{\cos x}{\left ( x^2+a^{2} \right )\left ( x^2+b^{2} \right )}~\mathrm{d}x$

First, I think you have a typo and $x$ should be $x^2$. In this case we can write $\displaystyle \frac{1}{(x^2+a^2)(x^2+b^2)}=\frac{1}{b^2-a^2}\left(\frac{1}{x^2+a^2}-\frac{1}{x^2+b^2}\right)$. Now you can split into two integrals of the form $\displaystyle \int\limits_{0}^{\infty}\frac{\cos(x)}{x^2+a^2}\text{d}x$ which you know how to compute.

For the general case of $\displaystyle \int\limits_{0}^{\infty}\frac{\cos(x)}{\displaystyle \prod_{i=1}^{\infty}(x^2+a_i^2)}\text{d}x$ I believe we can use the same method of breaking the multiplication into fractions of the form $\displaystyle \frac{1}{x^2+a_i^2}$ and using the known result.

Good luck.

Take $\;f(z):=\frac{e^{iz}}{z^2+a^2}\;,\;\;a>0\;$, and the contour

$$C_R:=[-R,R]\cup\Gamma_R\;,\;\;\Gamma_R:=\{z=Re^{it}\in\Bbb C\;;\;\;0<t<\pi\}\;,\;\;\;R\in\Bbb R^+$$

Observe that within the domain enclosed by the above contour our function has one unique simple pole: $\;z=ai\;$ , so

$$\oint_{C_R}f(z)\,dz=2\pi i\cdot\text{Res}(f)|_{z=ai}$$

Now

$$\text{Res}(f)|_{z=ai}=\lim_{z\to ai}(z-ai)f(z)=\frac{e^{-a}}{2ai}\implies$$

$$\frac\pi{ae^a}=2\pi i\frac{e^{-a}}{2ai}=\oint_{C_R}f(z)\,dz=\int_{-R}^R\frac{e^{ix}}{x^2+a^2}dx+\int_{\Gamma_R}f(z)\,dz$$

But

$$\left|\int_{\Gamma_R}f(z)\right|\le\ell(\Gamma_R)\cdot\max_{z\in\Gamma_R}\frac{e^{iz}}{z^2+a^2}\le\frac{\pi R}{2|(R^2e^{it}+a^2)|}\le\frac{\pi R}{2R^2-2a^2}\xrightarrow[R\to\infty]{}0$$

So finally

$$\frac\pi{ae^a}=\lim_{R\to\infty}\oint_{C_R}f(z)\,dz=\int_{-\infty}^\infty\frac{\cos x+i\sin x}{x^2+a^2}dx$$

Now just compare real parts (and observe the integrand is an even function).