With radar sounders, coherent backscattering simulations from global planetary DEMs typically display a deficit in diffuse clutter, which is mainly due to the implicit assumption that roughness at scales below the resolution of the DEM is absent. Indeed, while polynomial approximations of the phase evolution across the facet allow for fast and mathematically rigorous simulators, the coarse resolution of these planetary DEMs leads to a potentially significant portion of the backscattering response being neglected. In this paper, we derive the analytical phase response of a rough rectangular facet characterised by Gaussian roughness and a Gaussian isotropic correlation function under the linear phase approximation. Formulae for the coherent and incoherent power scattered by such an object are obtained for arbitrary bistatic scattering angles. Validation is done both in isolation and after inclusion in different Stratton-Chu simulators. In order to illustrate the different uses of such a formulation, we reproduce two lunar radargrams acquired by the LRS instrument with a Stratton-Chu simulator incorporating the proposed rough facet phase integral, and we show that the original radargrams are significantly better-reproduced than with state-of-the-art methods, at a similar computational cost. We also show how the rough facet integral formulation can be used in isolation to better characterise subglacial water bodies on Earth.