The goal of this work is to identify a density function of a physical body from a given data as the result of X-rays traveling through the body under different angles and offsets. The mathematical relation between parameter and data is described by the Radon transform. More specifically, the aim of the reconstruction is to find the singularity set of the density --- the collection of lines across which the density is discontinuous --- and the smooth density distribution on the complement of the singularity set.Mumford-Shah models are designed to extract simultaneously functional and geometric information for inaccessible parameters from indirect measurements. We propose a piecewise smooth Mumford-Shah model for the inversion of the tomography data. In our approach the functional variable is eliminated by solving a classical variational problem for each fixed geometry. The solution is then inserted in the Mumford-Shah cost functional leading to a geometrical optimization problem for the singularity set. The resulting shape optimization problem is solved using a shape sensitivity calculus and a propagation of shape variables in the level-set form.A specific difficulty poses the solution of the optimality system for the fixed geometry, which has the form of a coupled system of integro-differential equations on variable and irregular domains.A finite difference method based approach for the determination of a piecewise smooth density function as thesolution of the optimality system is presented. A second order accurate finite difference discretization is proposed.Here a standard five-point stencil is used on regular points of the underlying uniform grid and modificationsof the standard stencil are made at points close to the contour. The optimality system is solved iteratively.The Fast Multipole Method is used to solve domain integrals and numerical experiments are presented.