The Hermite-Gauss functions appear commonly in physics. These functions are formed from the product of a Hermite polynomial and a Gaussian:
$ u_n(x) = \left(\frac{2}{\pi w_0^2}\right)^{1/4} \frac{1}{\sqrt{n! 2^n}} H_n\left(\frac{\sqrt{2}x}{w_0}\right)\exp\left\{-\left(\frac{x}{w_0}\right)^2\right\}$
and are orthonormal:
$ \int_{-\infty}^{\infty} u_n(x) u_m(x) dx = \delta_{n,m}$
In a paper (2004 J. Opt. B 6 495) I found the following identity, which gives the decomposition of a displaced mode $u_0(x-a)$ in terms of a series over high-order Hermite-Gauss functions $u_n(x)$:
$ \int_{-\infty}^{\infty} u_0(x - a) u_n(x) dx = \frac{a^n}{w_0^n \sqrt{n!}} \exp\left\{ -\frac{a^2}{2 w_0^2}\right\} $
How is this derived?
(In addition to deriving it by hand, I would like to know how to coax Mathematica into giving it.)
EDIT: Here is an animation showing the decomposition of a displaced Gaussian into higher-order Hermite-Gauss functions (modes):