I want to calculate $V(\vec{x})=\int_{\mathbb{R}^3}\frac{\theta(\vec{y})}{|\vec{x}-\vec{y}|}dy$, where $\theta(\vec{y})=|\vec{y}|^2-R^2$ for $|\vec{y}|\leq R$ and $\theta(\vec{y})=0$ otherwise.
I can't get any further than recognising that $V(\vec{x})=\int_{|\vec{y}|\leq R}\frac{|\vec{y}|^2-R^2}{|\vec{x}-\vec{y}|}dy$. Switching to spherical co-ordinates makes sense but then I get a messy expression for $|\vec{x}-\vec{y}|$.
You have to choose the system of coordinates (directions of coordinate axes) while integrating in $\mathbb{R}^3$ - initially the Y system is not fixed in the space. You may choose the z axe in the direction of $\vec{x}$, and $|\vec{x}-\vec{y}|$ will be expressed as a function of the dot products of vectors, depending only on the polar ($\theta$) angle between $\vec{x}$ and $\vec{y}$ and not depending on radial angle $\varphi$. Integration over $\varphi$ is straightforward.
Your integral will look like $V(\vec{x})=2\pi\int_0^R \int_0^\pi\frac{y^2-R^2}{\sqrt{y^2+x^2-2xy\cos\theta}}y^2dy\sin\theta d\theta$, where $y$ and $x$ are $|\vec{y}|$ and $|\vec{x}|$ correspondingly. It is easier to integrate over $\theta$ first.
Hopefully it would help.