  Gravitational waves (GWs) influence the arrival times of radio signals coming
from pulsars. Here, we investigate the harmonic space approach to describing a
pulsar's response to GWs. We derive and discuss the "diagonalized form" of the
response, which is a sum of spin-2-weighted spherical harmonics of the GW
direction multiplied by normal (spin-weight 0) spherical harmonics of the
pulsar direction. We show how this allows many useful objects, for example, the
Hellings and Downs two-point function, to be easily calculated. The approach
also provides a clear description of the gauge dependence. We then employ this
harmonic approach to model the effects of angular correlations in the sky
locations of GW sources (sometimes called "statistical isotropy"). To do this,
we construct rotationally invariant ensembles made up of many Gaussian
subensembles, each of which breaks rotational invariance. Using harmonic
techniques, we compute the cosmic covariance and the total covariance of the
Hellings and Downs correlation in these models. The results may be used to
assess the impact of angular source correlations on the Hellings and Downs
correlation, and for optimal reconstruction of the Hellings and Downs curve in
models where GW sources have correlated sky locations.
