We present a principled study on establishing Gaussian processes over variables on the product of directional manifolds. As a basic functional component, a manifold-adaptive kernel is presented based on the von Mises distribution for Gaussian process regression on unit circles. Afterward, a novel hypertoroidal von Mises kernel is introduced to enable topology-aware Gaussian processes on hypertori with consideration of correlational circular components. Based thereon, we enable multi-output regression for learning vector-valued functions on hypertori using intrinsic coregionalization model and provide analytical derivatives in hyperparameter optimization. The proposed multi-output hypertoroidal Gaussian process is further embedded to a data-driven recursive estimation scheme for learning unknown range sensing models of angle-of-arrival inputs. Evaluations on range-based localization show that the proposed scheme enables superior tracking accuracy over parametric modeling and common Gaussian processes.