Earthquake ground motion spatial variability can influence significantly the response of certain structures. In order to accurately evaluate probabilistic characteristics of the seismic response of structures, the Monte Carlo simulation technique is still the only universal method of analysis when strong nonlinearities and input uncertainties are involved. Consequently, realizations of ground motion time histories taking into account both time and spatial variability need to be generated. Furthermore, for some design applications, the generated time histories must also satisfy the provision imposed by certain seismic codes stating that they have to be also response-spectrum-compatible. For these purposes, a spectral-representation-based methodology for generating fully non-stationary and spectrum-compatible ground motion vector processes at a number of locations on the ground surface is proposed in this paper. The simulated time histories do not require any iterations on the individual generated sample functions so that Gaussianity and prescribed coherence are suitably preserved. The methodology has also the advantage of providing the fully non-stationary and spectrum-compatible cross-spectral density matrix of the ground motion time-histories that can be used for reliability studies in an analytic stochastic fashion.