We consider the computation of the harmonic 3D spatial-domain Green's function of the surface of a piezoelectric substrate. Starting from the Green's function expressed in the spectral domain, the singular contributions are isolated and treated separately. It is found that the surface acoustic wave contributions give rise to an anisotropic generalization of the Hankel function H0(2), the spatial Green's function for the scalar two-dimensional wave equation. The asymptotic behavior at infinity and at the origin are also explicitely treated. The remaining non-singular part of the spatial Green's function is obtained numerically by a combination of fast Fourier transform and quadrature. Illustrations are given in the case of a substrate of Y + 128-cut lithium niobate.