We consider the inverse problem of determining spatially heterogeneous absorption and diffusion coefficients µ(x), D(x), from a single measurement of the absorbed energy E(x) = µ(x)u(x), where u satisfies the elliptic partial differential equationThis problem, which is central in quantitative photoacoustic tomography, is in general ill-posed since it admits an infinite number of solution pairs. Using similar ideas as in [31], we show that when the coefficients µ, D are known to be piecewise constant functions, a unique solution can be obtained. For the numerical determination of µ, D, we suggest a variational method based on an Ambrosio-Tortorelli approximation of a MumfordShah-like functional, which we implement numerically and test on simulated two-dimensional data.