How can I write a Gaussian state as a squeezed, displaced thermal state? I would like to write a Gaussian state with density matrix $\rho$ (single mode) as a squeezed, displaced thermal state:
\begin{gather}
\rho = \hat{S}(\zeta) \hat{D}(\alpha) \rho_{\bar{n}} \hat{D}^\dagger(\alpha) \hat{S}^\dagger(\zeta) .
\end{gather}
Here,
\begin{gather}
\rho_{\bar{n}} = \int_{\mathbb{C}} P_{\bar{n}}(\alpha) |\alpha\rangle\langle\alpha|d\alpha \text{ with } P_{\bar{n}}(\alpha) = \frac{1}{\pi \bar{n}} e^{- |\alpha|^2 / \bar{n}}
\end{gather}
is a thermal state with occupation $\bar{n}$,
\begin{gather}
\hat{S}(\zeta) = e^{(\zeta^* \hat{a}^2 + \zeta \hat{a}^{\dagger 2}) / 2}
\end{gather}
is the squeezing operator, and
\begin{gather}
\hat{D}(\alpha) = e^{\alpha^* \hat{a} - \alpha \hat{a}^\dagger} .
\end{gather}
is the displacement operator. I prefer to use the convention $\hat{x} = (\hat{a} + \hat{a}^\dagger) / \sqrt{2}$ and $\hat{p} = (\hat{a} - \hat{a}^\dagger) / \sqrt{2} i$.
I assume that the way to accomplish this is to derive the mean and variance of our Gaussian state $\rho$ and thereby determine $\zeta$ and $\alpha$. However, I have been unsuccessful in doing so. That is, given the mean and variance of our Gaussian state $\rho$, what are $\zeta$ and $\alpha$?
On a side note, I was also wondering if there is a standard result for the commutator of $\hat{S}(\zeta)$ and $\hat{D}(\alpha)$?
 A: I will follow the notes by A. Ferraro et al.
A state $\rho$ of a system with $n$ degrees of freedom is said to be Gaussian if its Wigner function can be written as
$$W[\rho](\boldsymbol{\alpha}) 
=
 \frac{\exp\left( -\frac{1}{2}(\boldsymbol{\alpha} - \bar{\boldsymbol{\alpha}})^T
\boldsymbol{\sigma}^{-1}_\alpha 
(\boldsymbol{\alpha} - \bar{\boldsymbol{\alpha}})
\right)
}{(2\pi)^n\sqrt{\text{Det}[\boldsymbol{\sigma}_\alpha ]}},$$
where $\boldsymbol{\alpha}$ and $\bar{\boldsymbol{\alpha}}$ are vectors containing all the $2n$ quadratures of the system and their average values, respectively, and $\boldsymbol{\alpha}$ is the covariance matrix, whose elements are defined as
$$[\boldsymbol{\sigma}]_{kl}
:= \frac{1}{2}
\langle \{R_k,R_l \} \rangle
-
\langle R_k\rangle  \langle R_l\rangle,
$$
where $\{\cdot,\cdot\}$ is the anticommutator, and $R_k$ is the $k-th$ element of the vector $\boldsymbol{R}= (q_1,p_1,\ldots,q_n,p_n)^T$ with the $q$s and $p$s being the position and momentum-like operators.
A very important result is that it turns out that Gaussian states can be fully characterised by their covariace matrix plus the vector of expectation values of the quadratures, $\bar{\boldsymbol{\alpha}}$. If your system only has one mode (one boson), then you only need a symmetric $2\times2$ matrix and two real numbers ($q$ and $p$) to describe it! This means a total of five parameters.
As you point out, we can write any Gaussian state as
$$
\rho = D(\bar{\alpha})S(\xi)\rho_{th}S^\dagger(\xi)D^\dagger(\bar{\alpha})
$$
where here $\bar{\alpha}\equiv\frac{1}{\sqrt{2}}(\bar{x}+i\bar{p})$ and $\xi = r e^{i\varphi}$. If your thermal state has mean photon number $N$, then it suffices to know $N$, $r$ and $\varphi$ to compute the covariance matrix. Its elements are given by
$$
\sigma_{11} =\frac{2N+1}{2} \left(\cosh (2r)+\sinh (2r)\cos (\varphi)\right)
$$
$$\sigma_{22} =\frac{2N+1}{2}\left( \cosh(2r)-\sinh(2r)\cos(\varphi) \right)$$
$$\sigma_{12} =\sigma_{21}=-\frac{2N+1}{2}\sinh(2r)\sin(\varphi).$$
You can see that the main properties of the state are captured by the covariance matrix, because the displacement $\bar{\alpha}$ can always be disregarded by local operations (it is a phase-space translation). In other words, you can always put it to zero. 
Answering your question, note that a Gaussian state is not simply a Gaussian distribution. You need more parameters than simply the variance and the mean (as you would to define a classical Gaussian probability distribution). These are in general five real values, but the essential ones are the ones entering the covariance matrix, as explained before.
As for the commutator, I am not aware of any closed formula. But I do know that displacing and then squeezing produces a state which has the same squeezing as the zqueezed-displaced state, but a different displacement.
