With each C-2-domain and material parameter, an infinite number of tensors, called the Generalized Polarization Tensors (GPTs), is associated. The GPTs contain significant information on the shape of the domain and its material parameter. They generalize the concept of Polarization Tensor (PT), which can be seen as the first-order GPT. It is known that given an arbitrary shape, one can find an equivalent ellipse or ellipsoid with the same PT. In this paper we consider the problem of recovering finer details of the shape of a given domain using higher-order polarization tensors. We design an optimization approach which solves the problem by minimizing a weighted discrepancy functional. In order to compute the shape derivative of this functional, we rigorously derive an asymptotic expansion of the perturbations of the GPTs that are due to a small deformation of the boundary of the domain. Our derivations are based on the theory of layer potentials. We perform some numerical experiments to demonstrate the validity and the limitations of the proposed method. The results clearly show that our approach is very promising in recovering fine shape details.