Let $\mu$ be a probability measure on $(\mathbb R^d,\mathcal B(\mathbb R^d))$, and let $\xi_{1/k}$ denote the $d$-fold product of the Gaussian distribution with mean zero and variance $\frac{1}{k}$. Let $\varphi$ denote the characteristic function of $\mu$. Am trying to show that the convolution $\mu \ast \xi_{1/k}$ has a density $f$ with respect to the Lebesgue measure $\lambda^d$, with
$$f(z)=\frac{1}{(2\pi)^d}\int_{\mathbb R^d} \varphi(x) e^{-i\langle z,x \rangle}e^{-\frac{\|x\|^2}{2k}}\lambda^d(dx).$$
I did the following:
Let $A\in \mathcal B(\mathbb R^d)$ and let $g_{1/k}$ denote the density of $\xi_{1/k}$. We have
$$(\mu \ast \xi_{1/k})(A)=\int_{\mathbb R^d\times\mathbb R^d} 1_{\{x+y\in A\}}d(\mu \otimes \xi_{1/k})\overset{Tonelli}{=}\int_{\mathbb R^d}\int_{\mathbb R^d} 1_{\{x+y\in A\}}\xi_{1/k}(dy)\mu(dx)=\int_{\mathbb R^d}\int_{\mathbb R^d} 1_{\{x+y\in A\}}g_{1/k}(y)\lambda^d(dy)\mu(dx)$$
$$\overset{z=x+y}{=}\int_{\mathbb R^d}\int_{\mathbb R^d} 1_{\{z\in A\}}g_{1/k}(z-x)\lambda^d(dz)\mu(dx)\overset{Tonelli}{=}\int_{\mathbb R^d}\int_{\mathbb R^d} 1_{\{z\in A\}}g_{1/k}(z-x)\mu(dx)\lambda^d(dz)=\int_A \int_{\mathbb R^d} g_{1/k}(z-x)\mu(dx)\lambda^d(dz),$$
which shows that $f(z)=\int_{\mathbb R^d} g_{1/k}(z-x)\mu(dx)$ is a density for $\mu \ast \xi_{1/k}$ with respect to $\lambda^d$. Since $\xi_{1/k}$ is a $d$-fold product of $N(0,\frac{1}{k})$, we have
$$ g_{1/k}(y)=\bigg(\frac{k}{2\pi}\bigg)^{d/2}\exp \bigg\{{\frac{-k\|y\|^2}{2}}\bigg\}.$$
Now, if $\phi$ denotes the characteristic function of $\xi_{k}$, then by using Fubini's theorem (complex case) and the known formula for the characteristic function of $N(0,k)$ we get
$$\phi(t)=\int_{\mathbb R^d}e^{i\langle t,y \rangle}\xi_k(dy)=\exp \bigg\{{\frac{-k\|t\|^2}{2}}\bigg\}.$$
Therefore, regarding $f(z)$ as a complex integral, we get
$$f(z)=\bigg(\frac{k}{2\pi}\bigg)^{d/2}\int_{\mathbb R^d} \int_{\mathbb R^d}e^{i\langle x-z,y \rangle}\xi_k(dy)\mu(dx).$$
Applying Fubini's theorem (complex case) once more, we obtain
$$f(z)=\bigg(\frac{k}{2\pi}\bigg)^{d/2}\int_{\mathbb R^d} \int_{\mathbb R^d}e^{i\langle x-z,y \rangle}\mu(dx)\xi_k(dy)=\bigg(\frac{k}{2\pi}\bigg)^{d/2}\int_{\mathbb R^d} \bigg[\int_{\mathbb R^d} e^{i\langle x,y \rangle}\mu(dx)\bigg] e^{-i\langle z,y \rangle}\xi_k(dy)$$
$$=\bigg(\frac{k}{2\pi}\bigg)^{d/2}\int_{\mathbb R^d} \varphi(y)e^{-i\langle z,y \rangle}\xi_k(dy)=\frac{1}{(2\pi)^d}\int_{\mathbb R^d} \varphi(y)e^{-i\langle z,y \rangle}\exp \bigg\{{\frac{-\|y\|^2}{2k}}\bigg\}\lambda^d(dy).$$
Is this correct? Thanks a lot for your comments/feedback.
This is a comment, that I think fits better here:
In general, suppose $\mu$, $\nu$ are finite Borel measures on $\mathbb{R}^n$ and that $\nu\ll\lambda_n$ ($\lambda_n$ is Lebesgue's measure) and that $g=\frac{d\nu}{d\lambda_n}$. Then $\int e^{-2\pi ix\cdot t}(\mu*\nu)(dx)=\widehat{\mu*\nu}(t)=\widehat{\mu}(t)\widehat{g}(t)$, where $\widehat{g}(t)=\int e^{-2\pi ix\cdot t}g(x)\,dx$ and $\widehat{\mu}(t)=\int e^{-2\pi ix\cdot t}\mu(dx)$.
From this equation, you can substitute $g$ for the density off a Gaussian (or any other density) to derive several identities. On the other hand $\mu*\nu\ll\lambda_n$ and $$f:=\frac{d\mu*\nu}{d\lambda_k}(x)=\int g(x-y)\mu(dy)$$ Hence $\widehat{f}=\widehat{\mu}\widehat{g}$.
If $\widehat{g}\in L_1$, then by Fourier's inversion formula $ g(x)=\widehat{\widehat{g}}(-x)$ and so $$f(y)=\int \Big(\int e^{2\pi iz\cdot(x-y)}\widehat{g}(z)\,dz\Big)\,\mu(dx)=\int e^{-2\pi iy\cdot z}\widehat{g}(z)\Big(\int e^{2\pi ix\cdot z}\mu(dx)\Big)\,dz$$