This paper deals with bistatic subsurface probing of a horizontally layered dielectric half-space by means of ultra-wideband electromagnetic pulses. A receiver collects reflections from the air-ground interface and from the gradients of dielectric permittivity in the half-space. This scenario is of interest for ground penetrating radar (GPR) applications. For the analytical description of the received signal, we developed and implemented a novel time-domain version of the coupled-wave Wentzel–Kramers–Brillouin approximation. Our solution is in very good agreement with finite-difference time-domain results, radically accelerates calculations, and effectively accounts for the protracted return signals observed in the lower part of the GPR spectrum. The paper includes results showing the application of the proposed technique to two case studies: in particular, the method was employed for the post-processing of experimental radargrams collected on Lake Chebarkul, in Russia, and to simulate GPR probing of the Moon surface, to detect smooth gradients of the dielectric permittivity in lunar regolith.