I have the following complex integral that corresponds to a complex integral of a Wigner function of a 1-mode Gaussian state:
$$I_n(\sigma) = \int^\infty_{-\infty} d^2\alpha \; L_n \left(\frac{4|\alpha|^2}{2 - \sigma}\right) \exp\left\{\frac{2\sigma |\alpha|^2}{\sigma - 2}\right\} \exp\left\{-\frac{1}{2}(\alpha - \lambda)^\top V^{-1}(\alpha - \lambda)\right\},$$
where $\alpha$ is a complex variable, $d^2\alpha = d\Re(\alpha)d\Im(\alpha)$, $\lambda$ and $V$ are the first (mean) and second (covariance) moments of the Gaussian state defined through:
$$\lambda = (a, b)^\top, \quad V = \begin{pmatrix} c & d \\ d & e \end{pmatrix},$$
where the coefficients $a,b,c,d,e$ are general complex numbers, and $L_n(x)$ are the Laguerre polynomials.
I am aware that this type of integral is usually evaluated by using polar coordinates. However, I am unsure how to achieve this.
Any help is appreciated. Even a method that obtains the solution numerically is fine.