Resolving multiple supermassive black hole binaries with pulsar timing arrays

We study the capability of a pulsar timing array (PTA) to individually resolve and localize in the sky monochromatic gravitational wave (GW) sources. Given a cosmological population of inspiralling massive black hole binaries, their observable signal in the PTA domain is expected to be a superposition of several nearly-monochromatic GWs of different strength. In each frequency bin, the signal is neither a stochastic background nor perfectly resolvable in its individual components. In this context, it is crucial to explore how the information encoded in the spatial distribution of the array of pulsars might help recovering the origin of the GW signal, by resolving individually and locating in the sky the strongest sources. In this paper we develop a maximum-likelihood based method finalized to this purpose. We test the algorithm against noiseless data showing that up to P/3 sources can be resolved and localized in the sky by an array of P pulsars. We validate the code by performing blind searches on both noiseless and noisy data containing an unknown number of signals with different strengths. Even without employing any proper search algorithm, our analysis procedure performs well, recovering and correctly locating in the sky almost all the injected sources.