I've been trying to prove by hand the Peskin's formula for the retarded propagator of the Klein Gordon equation, that is,
$$\int_{x^0 > y^0} \frac{d^4p}{(2\pi)^4} \frac{-e^{-ip(x-y)}}{i(p^2 - m^2)} = \int \frac{d^3p}{(2\pi)^3} \frac{e^{-ip(x-y)} - e^{ip(x-y)}}{2E_p}$$
With the condition that in the second integral, $p^0 = E_p$.
So I performed the standard contour integral of $\int \frac{dp^0}{2\pi} \frac{-e^{-ip(x-y)}}{i(p^2 - m^2)}$, contour that runs along the real axis except for the pole (where the contour is a small half circle around the pole), and then half circle back to the original position.
The integral for the half circle back seems indeed to go to $0$ in the infinite limit :
$$\int_0^\pi \frac{d\theta}{2\pi} (ire^{i\theta}) \frac{-e^{-i[re^{i\theta}(x^0-y^0) - \vec p(\vec x-\vec y)]}}{i((r e^{i\theta})^2 - \vec p^2 - m^2)} = \int_0^\pi \frac{d\theta}{2\pi} \frac{-e^{-i[re^{i\theta}(x^0-y^0) - \vec p(\vec x-\vec y)]}}{r e^{i\theta} - \frac{\vec p^2 + m^2}{r e^{i\theta}}}$$
Which in the limit of $r \rightarrow \infty$ will give something like $$\frac{-e^{-i[re^{i\theta}(x^0-y^0) - \vec p(\vec x-\vec y)]}}{r e^{i\theta} - 0}$$
And as the exponential is bounded, should go to 0.
When I computed around the poles, I used the following lemma :
$$\lim_{\epsilon\to0} \int_{\gamma_{\epsilon}} f(\zeta) \;d\zeta=i\theta_0 \mathrm{Res} \{ f(z), z_0 \}$$
The two poles are $p^0 = \pm E_p$, with residues $$ \frac{-1}{2\pi i}\frac{e^{-i[(\pm E_p)(x^0-y^0) - \vec p(\vec x-\vec y)]}}{\pm 2E_p}$$
So those two integrals are
$$I_\pm = \mp \frac{1}{4E_p}e^{-i[(\pm E_p)(x^0-y^0) - \vec p(\vec x-\vec y)]}$$
Meaning that the principal values of the integral over the real line should be
$$\text{P.V.}\int \frac{dp^0}{2\pi} \frac{-e^{-ip(x-y)}}{i(p^2 - m^2)} + I_+ + I_- + 0 = 0$$
$$\text{P.V.}\int \frac{dp^0}{2\pi} \frac{-e^{-ip(x-y)}}{i(p^2 - m^2)} = \frac{1}{4E_p}(e^{-i[E_p(x^0-y^0) - \vec p(\vec x-\vec y)]} - e^{i[ E_p(x^0-y^0)+ \vec p(\vec x-\vec y)]})$$
Which has a factor of half in excess compared to the original (those things happen tho), but more importantly contains $E_p(x^0-y^0)+ \vec p(\vec x-\vec y)$ in its second exponential instead of $E_p(x^0-y^0)- \vec p(\vec x-\vec y)$. Not quite sure where that change should come from, as nowhere did I touch the spatial part of the momentum in the calculations.
