In the case of a large class of static, spherically symmetric black hole solutions in higher order modified gravity models, an expression for the associated energy is proposed and identified with a quantity proportional to the constant of integration, which appears in the explicit solution. The identification is achieved making use of derivation of the First Law of black hole thermodynamics from the equations of motion, evaluating independently the entropy via Wald method and the Hawking temperature via quantum mechanical methods in curved space-times. Several non trivial examples are discussed, including a new topological higher derivative black hole solution, and the proposal is shown to work in all examples considered.