Evaluation of an integral of the form $C({\bf r})=\int\frac{{\rm d}^dk}{(2\pi)^d}\frac{{\rm e}^{{\rm i}{\bf k}\cdot{\bf r}}}{|{\bf k}|^2+m^2}$ Consider a momentum-integral of the form in $d$ spatial dimensions:$$C({\bf r})=\int\frac{{\rm d}^dk}{(2\pi)^d}\frac{{\rm e}^{{\rm i}{\bf k}\cdot{\bf r}}}{|{\bf k}|^2+m^2}$$where ${\bf r}=(x_1,x_2,\cdots,x_d)$, ${\bf k}=(k_1,k_2,\cdots,k_d)$ and $m^2$ is a real number, could be positive, negative or zero. (Despite my notation, $m^2$ is any constant and has nothing to do with mass) See Eq.$(4.9)$ of this lecture note.
Question Is this integral doable? I am interested in an exact analytical formula for extracting the results for different spatial dimensions $d$ (namely, $1,2$ and $3$).

Update 
The existing answer by @DinosaurEgg, tells how to compute $C(\bf r)$ for $m^2>0$. I would appreciate if someone can ameliorate expand on his answer by addressing how to handle the case $m^2<0$. 
 A: For $m^2>0$, the integral clearly converges, so this is the case I will tackle. Using the identity
$$\frac{1}{\mathbf{k}^2+m^2}=\int_{0}^{\infty}ds~e^{-s(\mathbf{k}^2+m^2)}$$
we can rewrite the integral as follows
$$\begin{align}C(\mathbf{r})&=\int_{0}^{\infty}ds~e^{-sm^2}~\prod_{n=1}^{d}\int_{-\infty}^{\infty}\frac{dk_n}{2\pi}e^{ik_nr_n}e^{-sk_n^2}\\&=\int_{0}^{\infty}ds~e^{-sm^2}\Big(\frac{1}{2\sqrt{\pi s}}\Big)^de^{-\mathbf{r}^2/4s}\\&=\frac{1}{(4\pi)^{d/2}}\int_{0}^{\infty}ds~\frac{e^{-sm^2}e^{-\mathbf{r}^2/4s}}{s^{d/2}}\end{align}$$
which in turn is a well known integral representation of a modified Bessel function:
$$C(\mathbf{r})=\frac{m^{d-2}}{(2\pi)^{d/2}}\frac{1}{(mr)^{d/2-1}}K_{d/2-1}(mr)$$
@joriki's answer covers the ambiguity encountered when $m^2<0$. When $m^2=0$ and $d\geq 3$ this integral converges and reduces to the propagator of the Laplace equation:
$$C_{m=0}(\mathbf{r})=-\frac{1}{A_d(d-2)|\mathbf{r}|^{d-2}}$$
where $A_d$ denotes the surface area of the $d$-dimensional sphere.
EDIT:
I feel that this bounty went a little wasted, so let me elaborate a bit on the integral with $m^2<0$. Clearly, the integral as presented is nonsensical, as it passes through a pole at $\mathbf{k}^2=m^2$. However it makes sense to ask the following:
Question: What is the Green's function of the following Helmholtz equation in n-dimensions
$$\nabla^2 G- \lambda G=-\delta(\mathbf{r}-\mathbf{r}')$$
If we assume that $G$ is well behaved at infinity we can try to do a Fourier transform, and we would recover the exact integral quoted in the OP. Clearly, this is an assumption, and indeed it fails to be true whenever $\lambda=-m^2<0$, because the solutions of the homogeneous equation  then are oscillatory exponential-like  and hence the choice of solution is no longer uniquely specified, just like it is when $\lambda=m^2>0$ by the requirement that the Green's function decays quickly enough at infinity so that the integral converges.
An easy example to demonstrate that is $n=1$, where we can construct the solution easily by using the solutions of the homogeneous equation piecewise around the point $r'$ imposing the following conditions:
$$\Big(\frac{d^2}{dr^2}+m^2\Big)G(r,r')=-\delta(r-r')$$
$$G(r,r')=\begin{Bmatrix}G_+ &r-r'>0\\G_-&r-r'<0\end{Bmatrix}$$
$$\lim_{r\to r'^+}G_+=\lim_{r\to r'^-}G_-~~,~~\lim_{r\to r'^{+}}\frac{dG_+}{dr}-\lim_{r\to r'^{-}}\frac{dG_-}{dr}=-1$$ and the most general solution is:
$$G(r,r')=\begin{Bmatrix}A\cos m(r-r')+(B-\frac{1}{m})\sin m(r-r') &r-r'>0\\A\cos m(r-r')+B\sin m(r-r')&r-r'<0\end{Bmatrix}$$
This has a two parameter family of solutions, all of which have a different Fourier transform. If one has to think of the integral above, we can only think of it as an object that represents a class of Green's functions, and then we would have to find a way to extract them from it. See a reference here for a methodology for extracting propagators with specific asymptotic properties from the integral.
For example, we can obtain one of the Green's functions closely related to the Feynman propagator by giving the mass a small negative imaginary part. Computing the integral using residue calculus we get
$$\lim_{\epsilon\to 0^+}\int_{-\infty}^{\infty}\frac{dk}{2\pi}\frac{e^{ik(r-r')}}{k^2-m^2-i\epsilon}=i\Big(\frac{e^{im(r-r')}}{2m}\theta(r-r')+\frac{e^{im(-(r-r'))}}{2m}\theta(-(r-r'))\Big)=i\frac{e^{im|r-r'|}}{2m}$$
which is indeed the only solution that behaves like an outgoing wave and one can see that from the general Green's formula.
In higher dimensions the situation is not much different. For brevity we will assume that $G(\mathbf{r},\mathbf{r}')=u(|\mathbf{r}-\mathbf{r}'|)$ which is correct and can be shown and knowing that this is the functional form of the Green's function, we can show that it satisfies the singular delta equation
$$\mathcal{L} u=\frac{1}{r^{n-1}}\frac{d}{dr}\Big(r^{n-1}\frac{du}{dr}\Big)+m^2 u=-\frac{\delta(r)}{r^{n-1}A_n}$$
where $A_n$ is the surface area of the $n$-dimensional sphere. It can be shown that the equation can be transformed into the following inhomogeneous Bessel equation for $f=r^{n/2-1}u$
$$r^2\frac{d^2f}{dr^2}+r\frac{df}{dr}+(m^2r^2-(\frac{n}{2}-1)^2)f=-\frac{\delta(r)}{r^{n/2-2}A_n}$$
This is an inhomogeneous Bessel equation and thus the general solution for $G$ should be of the form
$$G(\mathbf{r},\mathbf{r}')=C_1\frac{J_{n/2-1}(m|\mathbf{r}-\mathbf{r}'|)}{|\mathbf{r}-\mathbf{r}'|^{n/2-1}}+C_2\frac{Y_{n/2-1}(m|\mathbf{r}-\mathbf{r}'|)}{|\mathbf{r}-\mathbf{r}'|^{n/2-1}}$$
However the $C_2$ term is fixed. It can be shown that
$$(\nabla^{2}+m^2)\frac{J_{n/2-1}(mr)}{r^{n/2-1}}=0$$
$$(\nabla^{2}+m^2)\frac{Y_{n/2-1}(mr)}{r^{n/2-1}}=C\frac{\delta(r)}{r^{n-1}}$$
In fact we can show that if we integrate on a small ball of radius $\epsilon$ and take the limit $\epsilon \to 0$ and using the asymptotic expansion of the Neumann function as $r\to 0$, we have
$$\begin{align}\int_{r\leq \epsilon}(\nabla^2+m^2)\frac{Y_{n/2-1}(mr)}{r^{n/2-1}}d^n\mathbf{r}&=\int_{r=\epsilon}\nabla\Big(\frac{-2^{n/2-1}\Gamma(\frac{n}{2}-1)}{r^{n-2}}\Big)\cdot d\mathbf{S}+\mathcal{O}(\epsilon^2)\\&=-(n-2)A_n2^{n/2-1}\Gamma(n/2-1)+\mathcal{O}(\epsilon^2)\end{align}$$
which determines the constant to be
$$C=-(n-2)(2/m)^{n/2-1}\Gamma(n/2-1)~~~,~~~ C_2=-\frac{1}{CA_n}=\frac{m^{n/2-1}}{2(2\pi)^{n/2}}$$
Finally, ingoing and outgoing solutions may be constructed using the Hankel functions but as before, with a constraint:
$$G(r)=i(A~ H^{(1)}_{n/2-1}(mr)+B~ H^{(2)}_{n/2-1}(mr))$$
$$A-B=C_2$$
One obtains the ingoing/outgoing solutions setting $B=0$/$A=0$ respectively.
A: This integral is not well-defined because of the singularity at $k=m$. It’s a relativistic propagator for which you can get different results with different physical meanings (retarded, advanced or Feynman propagators) depending on how you deform the contour in the complex plane to avoid the singularity; see Wikipedia.
