My computations seem to be correct but I get a contradiction. Considering the Gaussian in $\mathbb{R}^n$ given by $G_n(x)=\frac{1}{(2\pi)^{n/2}} e^{-x^2/2}$, we can decompose its Fourier transform as:
$\hat{G}_n(\lambda)=\frac{1}{(2\pi)^n}\int_{\mathbb{R}^n} e^{-x^2/2-ix\lambda}dx = \frac{1}{(2\pi)^{n/2}}\prod_{j=1}^n \frac{1}{\sqrt{2\pi}}\int_{\mathbb{R}}e^{-x_j^2/2-ix_j\lambda_j}dx_j = \frac{1}{(2\pi)^{n/2}} \prod_j \hat{G}_1(\lambda_j).$
Then, I compute $\hat{G}_1$ and find $\hat{G}_1(\lambda)= \frac{1}{\sqrt{2\pi}} \int_{\mathbb{R}} e^{- \frac{(x+i\lambda)^2}{2}}e^{-\lambda^2/2}dx = G_1(\lambda) \int_{\mathbb{R}} e^{- \frac{(x+i\lambda)^2}{2}}dx = G_1(\lambda) \int_{\mathbb{R}} e^{- \frac{x^2}{2}}dx = \sqrt{2\pi}G_1(\lambda).$
Putting this result back in the first equations gives $\hat{G}_n(\lambda)=\frac{1}{(2\pi)^{n/2}} \prod_j \sqrt{2\pi}G_1(\lambda_j)= \frac{1}{(2\pi)^{n/2}}e^{-\lambda^2/2}=G_n(\lambda).$
To me, there must be something I'm doing wrong since setting n=1 leads to a contradiction. Some help would be appreciated.
Please note that you are using this convention of Fourier transform: $$\hat{f}(\lambda) = \frac{1}{(2\pi)^{n/2}} \int_{\mathbb{R}^n} f(x) \exp(-ix\lambda) dx$$ Under this convention, the standard n-dimensional Gaussian distribution is invariant under the transform. So, there's a mistake in your computation of $\hat{G}_1$. It should be $$\hat{G}_1(\lambda)= \frac{1}{\sqrt{2\pi}} \frac{1}{\sqrt{2\pi}} \int_{\mathbb{R}} \exp \left(\frac{-x^2}{2} - ix\lambda \right)dx = \frac{1}{2\pi} \int_{\mathbb{R}} e^{- \frac{(x+i\lambda)^2}{2}}e^{-\lambda^2/2}dx$$