Integrating a matrix function involving a determinant and exponential trace

As you pointed out, the integrand is invariant under the action of the orthogonal group. Thus it suffices to integrate over diagonal matrices, and multiply the result by the volume of the orthogonal group.

Suppose $X=diag(x_1,\ldots,x_k)$. Then the integral becomes $$ \int_{\mathbb R} \prod_{i=1}^k\left(1+\frac{g}{a}x_i^2\right)^{-d/2}e^{-ax_i^2/2}\ dx_1\cdots dx_k=\left[\int_{\mathbb R}\left(1+\frac{g}{a}x^2\right)^{-d/2}e^{-ax^2/2}\ dx\right]^k. $$ The innermost integral may be expressed in terms of moments of the Gaussian distribution.