The equation that you posted should be the same as the one used to calculate the fidelity in the Gaussian backend. The function used to calculate it can be found here, in a package called The Walrus. As stated in the docstring, it’s an implementation of the square of Eq. 96 in ‘Gaussian states and operations - a quick reference’, Brask.
Comparing Eq. 101 in Brask, for single mode Gaussian states, with Eq. 80 in ‘Gaussian Quantum Information’, Weedbrook, et al., we should see (using some substitutions) that they are identical.
Let me know if this answers your question!