The dynamic behavior of jointed assemblies exhibiting friction nonlinearities features amplitude-dependent dissipation and stiffness. To develop numerical simulations for predictive and design purposes, macro-scale High Fidelity Models (HFMs) of the contact interfaces are required. However, the high computational cost of such HFMs impedes the feasibility of the simulations. To this end, we propose a model-driven method for constructing hyper-reduced order models of such assemblies. Focusing on steady-state analysis, we use the Multi-Harmonic Balance Method (MHBM) to formulate the equations of motion in frequency domain. The reduction basis is constructed through solving a set of vibration problems corresponding to fictitious interface conditions. Subsequently, a Galerkin projection reduces the order of the model. Nonetheless, the necessary fine discretization of the interfaces represents a bottleneck for achieving high speedups. For this reason, we implement an adapted Energy Conserving Weighing and Sampling (ECSW) technique for Hyper Reduction (HR), thereby allowing significant speedups for meshes of arbitrary fineness. This feature is particularly advantageous since analysts typically encounter a trade-off between accuracy and computational cost when deciding on the mesh size, whose estimation is particularly challenging for problems of this type. To assess the accuracy of our method without resorting to the HF solution, we propose an error indicator with thresholds that have proven reliable in our analyses. Finally, the accuracy and efficiency of the method are demonstrated by two case studies.