Scenario: I have an integral originally expressed in Cartesian coordinates which should--in principle--converge, but I have difficulty evaluating it explicitly (due to apparent divergences) in spherical coordinates using the residue theorem.
Setup: I want to evaluate the following 4-dimensional integral in Euclidean space,
\begin{equation} \mathcal{I}_L(x_1,x_2)=\int_{|y-x_2|\le L} d^4y \frac{1}{(y-x_1)^2 (y-x_2)^2}, \tag{1}\label{1} \end{equation}
where $y,x_1,x_2$ are 4-vectors; "$||$" denotes the Euclidean distance function; and $(y-x_1)^2\equiv(y-x_1)\cdot(y-x_1)$ with "$\cdot$" the dot product. Also assume $L>0$ but finite.
First note that although the integrand diverges at $y=x_1,x_2$, the integral should be expected to converge on dimensional grounds.
Also note that shifting the integration variable, $y\to y+x_2$, in \eqref{1} gives,
\begin{align} \mathcal{I}_L(x_1,x_2)=\int_{|y|\le L} d^4y \frac{1}{(y-(x_1-x_2))^2 y^2} \tag{2}\label{2}\\ \end{align}
Now this 4-dim'l integral can be simplified to a 2-dim'l integral by switching to 4-dim'l "spherical coordinates",
\begin{align}
y^1&=|y|\cos(\phi_1)\\
y^2&=|y|\sin(\phi_1)\cos(\phi_2)\\
y^3&=|y|\sin(\phi_1)\sin(\phi_2)\cos(\phi_3)\\
y^4&=|y|\sin(\phi_1)\sin(\phi_2)\sin(\phi_3),\\
\end{align}
where $\phi_1,\phi_2$ range over $[0,\pi]$ while $\phi_3$ ranges over $[0,2\pi]$. [Notation: superscripts denote components of Cartesian 4-vector--i.e. $y=(y^1,y^2,y^3,y^4)$]
Choosing to orient our axes such that $\phi_3$ coincides with the angle between the 4-vectors $y$ and $x_1-x_2$ allows us to write \eqref{2} as
\begin{align} \mathcal{I}_L(\sigma)&= \int_0^L d|y|\,|y|^3 \int_{S^3} d\Omega_3 \frac{1}{(|y|^2-2|y|\sigma \cos{\phi_3}+\sigma^2)|y|^2}\\ &=\pi \int_0^L d|y| \int_0^{2\pi} d\phi_3 \frac{|y|}{|y|^2-2|y|\sigma \cos{\phi_3}+\sigma^2}, \label{3}\tag{3}\\ \end{align} where I define the shorthand $\sigma\equiv |x_1-x_2|$.
Now, the angular part of the integral in \eqref{3} looks like a good candidate for a calculus of residues approach, so I take $z=e^{i\phi}$ and find
\begin{align} \mathcal{I}_L(\sigma)&=\pi i \int_0^L d|y| \oint_{|z|=1}dz \frac{|y| }{(|y|z-\sigma)(\sigma z-|y|)}\\ &=\frac{\pi i}{\sigma}\int_0^L d|y|\oint_{|z|=1}dz\frac{1}{(z-\sigma/|y|)(z-|y|/\sigma)}.\label{4}\tag{4}\\ \end{align}
Using the residue theorem, I find
\begin{align} \oint_{|z|=1}dz\frac{dz}{(z-\sigma/|y|)(z-|y|/\sigma)}=2\pi i \frac{\sigma |y|}{|y|^2-\sigma^2}\left\{ \theta (\sigma-|y|)-\theta (|y|-\sigma)\right\}, \tag{5}\label{5} \end{align}
where $\theta(x>0)=1$ and $\theta(x)=0$ otherwise.
Plugging \eqref{5} into \eqref{4}, yields
\begin{align} \mathcal{I}_L(\sigma)=-2\pi^2\int_0^L d|y|\frac{|y|}{|y|^2-\sigma^2}\left\{ \theta (\sigma-|y|)-\theta (|y|-\sigma)\right\} \label{*} \tag{*} \end{align}
Confusion/Contradiction: For $L<\sigma$, \eqref{*} does seem to give a convergent result: \begin{align} \mathcal{I}_{L<\sigma}(\sigma)=-\pi^2 \ln\left[1-L^2/\sigma^2\right].\\ \end{align} [Recall: I have defined the shorthand $\sigma\equiv |x_1-x_2|$]
However, for $L\ge\sigma$, $\eqref{*}$ seems to diverge which contradicts the obvious convergence of $\eqref{1}$.
More explicitly, for $L\ge\sigma$, the integral \eqref{*} naturally splits into 2 parts:
\begin{align} \mathcal{I}_L(\sigma)&=-2\pi^2\left[ \int_0^\sigma-\int_\sigma^L\right] d|y|\frac{|y|}{|y|^2-\sigma^2}\\ &=-\pi^2 \left[ \ln{(|y|^2-\sigma^2)}\right]^\sigma_0+\pi^2\left[\ln{(|y|^2-\sigma^2)}\right]^L_\sigma\\ &\sim -\infty\\ \end{align}
Edit: (elaboration on H.H.Rugh's answer) As correctly pointed out by H.H.Rugh below, $\phi_3$ cannot be chosen to coincide with the angle between $y$ and $x_1-x_2$. The nature of my mistake can be easily visualized in 3-dimensional Euclidean space, where it is conventional to choose one angle $\phi:[0,2\pi]$ to denote rotations in the "$x,y$"-plane around the "$z$-axis" and one angle $\theta:[0,\pi]$ to denote rotations away from the $z$-axis as follows:
\begin{align} x&=r \sin\theta \cos\phi\\ y&=r \sin\theta \sin\phi\\ z&=r \cos\theta,\\ \end{align} where $r^2=x^2+y^2+z^2$.
Clearly the $\theta$-component of an arbitrary 3-vector in spherical coordinates, $(r,\phi,\theta)$ coincides with the angle between that 3-vector and one lying along the $z$-axis. However, the $\phi$-component of that same arbitrary 3-vector does not necessarily coincide with the angle between that 3-vector and one lying along the $x$ or $y$-axis. This is analogous to my mistake of choosing $\phi_3$ to coincide with the angle between the 4-vectors $y$ and $x_1-x_2$.
For completeness, I also explicitly finish the computation started by H.H.Rugh, beginning with
\begin{align} \mathcal{I}_L(\sigma)&= 2\pi \int_0^L dr \int_0^{2\pi} d\phi \frac{r \sin^2 \phi}{r^2 - 2 r \sigma \cos \phi_1 + \sigma^2 } \\ &=\frac{i\pi}{2\sigma}\int_0^L dr \oint_{z=1} dz \frac{(z^2-1)^2}{z^2(z-\sigma/r)(z-r/\sigma)}\\ \end{align}
There is now an additional relevant residue at $z=0$. Applying the residue theorem, I find
$$ \mathcal{I}_{L}(\sigma)= \begin{cases} \pi^2 L^2/\sigma^2 & \text{if } L<\sigma \\ \pi^2 \left(1+\ln(L^2/\sigma)\right) & \text{if } L\ge \sigma\\ \end{cases} $$
I think that the divergence comes from a wrong setup for the spherical coordinates relative to $\Delta x=x_1-x_2$. I don't see how you can make $\phi_3$ the angle between $y$ and $\Delta x$? Instead it seems natural to let $\Delta x$ be in the (positive) $y^1$ direction so the angle is $\phi_1$.
The volume element in 4-dim spherical coordinates is $$r^3 \sin^2(\phi_1) \sin(\phi_2) dr \; d\phi_1\; d\phi_2\; d\phi_3$$ and the square distance between $y$ and $\Delta x$ is (with $\sigma=|\Delta x|$ and the described choice of angle) $$ (r \cos \phi_1 - \sigma )^2 + (r \sin \phi_1)^2 = r^2 - 2 r \sigma \cos \phi_1 + \sigma^2 $$ For $r=\sigma$ this in the denominator indeed gives rise to a singularity as $\phi_1$ goes to zero of order $1/\phi_1^2$ but when you integrate you have now (due to dim 4) the important factor $(\sin \phi_1)^2$ which will remove the singularity. The integral after carrying out the trivial part should look something like (modulo errors in my calculations) $$ 4\pi \int_0^\pi dr \int_0^\pi d\phi \frac{r \sin^2 \phi}{(r-\sigma)^2+ 4r\sigma (1-\cos\phi)} $$
You may replace by half the integral over $[0,2\pi]$ and try to evaluate using residues but in any case the integrals are now integrals of bounded functions.