We revisit the problem of predicting the spin magnitude and direction of the black hole (BH) resulting from the merger of two BHs with arbitrary masses and spins inspiraling in quasi-circular orbits. We do this by analyzing a catalog of 619 recent numerical-relativity simulations collected from the literature and spanning a large variety of initial conditions. By combining information from the post-Newtonian approximation, the extreme mass-ratio limit, and perturbative calculations, we improve our previously proposed phenomenological formulae for the final remnant spin. In contrast with alternative suggestions in the literature, and in analogy with our previous expressions, the new formula is a simple algebraic function of the initial system parameters and is not restricted to binaries with spins aligned/anti-aligned with the orbital angular momentum but can be employed for fully generic binaries. The accuracy of the new expression is significantly improved, especially for almost extremal progenitor spins and for small mass ratios, yielding an rms error sigma approximate to 0.002 for aligned/anti-aligned binaries and sigma = 0.006 for generic binaries. Our new formula is suitable for cosmological applications and can be employed robustly in the analysis of the gravitational waveforms from advanced interferometric detectors.
