Thermodynamic properties of ZrC are calculated up to the melting point (T-melt approximate to 3700 K), using density functional theory (DFT) to obtain the fully anharmonic vibrational contribution, and including electronic excitations. A significant improvement is found in comparison to results calculated within the quasiharmonic approximation. The calculated thermal expansion is in better agreement with experiment and the heat capacity reproduces rather closely a CALPHAD estimate. The calculations are presented as an application of a development of the upsampled thermodynamic integration using Langevin dynamics (UP-TILD) approach. This development, referred to here as two-stage upsampled thermodynamic integration using Langevin dynamics (TU-TILD), is the inclusion of tailored interatomic potentials to characterize an intermediate reference state of anharmonic vibrations on a two-stage path of thermodynamic integration between the original DFT quasiharmonic free energy and the fully anharmonic DFT free energy. This approach greatly accelerates the convergence of the calculation, giving a factor of improvement in efficiency of similar to 50 in the present case compared to the original UP-TILD approach, and it can be applied to a wide range of materials.