I will present the HEPfit framework, a flexible open-source C++ code to calculate observables and to fit various models in high energy physics. The Two-Higgs-Doublet models with a softly broken Z2 symmetry belong to these models. I will explain the different constraints coming from experiment (including LHC run I) and theory. Concerning the latter, I will focus on next-to-leading order effects in renormalization group running and in the unitarity of scalar-scalar to scalar-scalar scattering. Eventually, I will show the results of a combined HEPfit analysis in these models.