We proposed an approach for the estimation of the depth profile of thermally inhomogeneous samples from combined laterally and frequency resolved photothermal measurements. The mathematical procedure of data inversion makes use of the quasianalytical solution of the forward problem. Introducing an appropriately chosen grid of depth coordinates the Hankel transform of the surface temperature can be expressed by a continuous fraction formula. This enables the usage of the effective conjugated gradient technique to retrieve the thermal depth profiles by minimization of the objective function. Making use of the a priori information about the inhomogeneous sample we chose an appropriate Tikhonov‘s stabilizer function and by this way remarkably improved the iteration procedure. The success of this method consisted both in a drastically reduced computation time and in a better approximation of the searched profiles. This was demonstrated by numerical simulations modeling the case of surface hardened steel.