Furnas Volcano (São Miguel Island, Azores) is a steep-sided caldera complex resulting from several collapses. In this paper, we analyse and interpret the gravity survey carried out by the Instituto de Astronomı́a y Geodesia (Spain) and the gravity data published by the Instituto Nacional de Meteorologia e Geofisica (Portugal). The resulting Bouguer anomaly reveals the presence of a regional trend, several local features and a component of noise. In separating the regional trend we use a robust polynomial fitting. Further, a covariance analysis and a least squares prediction are applied to filter the noise and to model the local anomaly. A first signal map represents the main local anomaly, while the secondary signal map shows superficial anomalies. This distribution reveals a local main minimum close to the 1630 A.D. and Gaspar craters. Several possible alignments are shown that lead to a model, mainly with SW–NE and SE–NW orientations, of the negative anomaly, with the pattern of the low density zones showing a relationship with the volcanic features. Applying a stabilized linear method of gravimetric inversion, we obtain a density contrast distribution of the sources of the filtered local Bouguer anomaly, considering as unknowns the anomalous densities of a fixed discrete distribution of rectangular prisms which represents the subsoil volume. A mixed minimization criteria for the weighted residuals and the weighted density contrast allow us to obtain a solution as a discrete three-dimensional density contrast model. The obtained anomalous densities model shows a clear relation with the structure of the caldera complex and tectonic trends, with the principal body of negative contrast density being associated with the subsurface structure of the caldera.