How do I compute the eigenfunctions of the Fourier Transform?

One method is to apply the Bargmann transform:

$$B_f(z) = \int_{-\infty}^{\infty} f(x) \exp \left(2 \pi x z - \pi x^2 - \frac{\pi}{2}z^2 \right) dy$$

to the eigenvalue equation:

$$f(x) = \lambda \int_{-\infty}^{\infty}f(y) \exp(-2 \pi ix y) dy$$

to obtain the following relation for the Bargmann transformed eigenfunctions:

$$B_f(z) = \lambda B_f(-iz) $$

whose solutions are the monomials:

$B_f(z) = z^n$, corresponding to the eigenvalues $i^{n \bmod{4}}$.

Now, the inverse Bargmann transform of the monomials are just the Hermite functions.

Remark: The application of the Bargmann transform to the eigenvalue equation entails a change of the integration order which is possible because the Hermite functions are bounded.

The projection operators onto the $n$th subspace in the Bargamann representation are given by the kernels

$$ P(v,\bar{z})=\frac{(v\bar{z})^n}{n!} ,$$

which act on the Bargmann transformed functions according to:

$$ Pf(v) = \int_{\mathbb{C}}P(v,\bar{z}) f(z) \exp(-z\bar{z})dzd\bar{z} .$$

In the time domain the projection kernels are just $P(t, \tau) = H_n(\tau) H_n(t)$ due to the orthonormality of the Hermite functions.