He aquí cómo Euler hizo esto (de E675, traducción disponible aquí). Empezamos con
$$ \int_0^{\infty} x^{m-1} e^{-x} \, dx = \Gamma(m). $$
El cambio de las variables de da
$$ \int_0^{\infty} x^{m-1} e^{-kx} \, dx = \frac{\Gamma(m)}{k^m}. $$
Euler ahora se supone que esto todavía funciona si $k=p \pm iq$ es complejo, siempre $p>0$. (Que no, pero esto necesita una aplicación de Cauchy teorema.) Entonces tenemos
$$ \int_0^{\infty} x^{m-1} e^{-(p \pm iq)x} \, dx = \frac{\Gamma(m)}{(p \pm iq)^m}, $$
y si escribimos $p=f\cos{\theta}$, $q=f\sin{\theta}$, podemos aplicar la fórmula de Euler para obtener
$$ \int_0^{\infty} x^{m-1} e^{-(p \pm iq)x} \, dx = \frac{\Gamma(m)}{f^m}(\cos{m\theta} \mp i\sin{m\theta}). $$
Tomando la parte real nos da
$$ \int_0^{\infty} x^{m-1} e^{-px} \cos{qx} \, dx = \frac{\Gamma(m)}{f^m}\cos{m\theta}, $$
y el resultado sigue invirtiendo las expresiones de $p$ $q$ en términos de$f$$\theta$.