Geophysical geometry inversion aims to reconstruct the geometrical characteristics of subsurface target bodies, which is different from conventional inversion techniques that focus on subsurface physical properties (e.g., density and velocity). The published works on geometry inversion either focus on recovering one “optimal” model without considering uncertainty or are limited to situations where only a single anomalous body is sought when uncertainty is quantified. We aim to recover the geometries of multiple anomalous bodies and quantify their uncertainty in the trans-dimensional Bayesian framework. We adopt a sparse geometry parameterization strategy which allows approximating the shapes of anomalous bodies via a set of simple geometries (e.g., ellipses). To address the problem of multiple bodies in an area of study, we propose two strategies. The first strategy randomly samples a tuning parameter associated with the computation of the alpha shapes. In the second strategy, we devise two categories of geometries for sampling parent and child geometries. Each child geometry is assigned to one parent geometry, and various geometrical families represent the multiple target bodies. To understand the effects of multiple data constraints and structural constraints, we design several scenarios, where the borehole gravity and gravity gradient data serve as additional data constraints and seismic imaging provides the structural constraint on the top of target bodies. We apply our method to the modified Pluto salt model consisting of three salt bodies with drastically different shapes and volumes. The results demonstrate that our framework is capable of recovering the disconnected anomalous bodies while taking into account multiple data constraints as well as structural constraints.