The gravity method is one of the geophysical methods that has been used in arnwider range of geophysical prospecting and investigations. One of the indispensablernsteps in this method is the inverse modelling of the measured data to estimate thernsubsurface density distribution and geometrical properties (e.g. shape and depth) of therncausative bodies. Particularly, delineating localized and blocky geologic features,rnthrough the inversion of gravity data is an important goal in a range of geophysica lrninvestigations and it is still a subject of interest and concern in the scientific communit y.rnMost conventional inversion algorithms generally yield smooth models with poor edgerndefinition. These methods have difficulties in recovering non-smooth distributions thatrnhave sharp boundaries. The main objective of this thesis was then to develop andrnimplement a gravity inversion algorithm that can produce compact and sharp images,rnaiming at recovering localized and blocky geologic features with varying geometricrnrepresentations. In the course of the thesis work, this goal has been achieved byrnpresenting a gravimetric inverse modelling method that has been tested to be effective.rnAt the heart of the developed inversion method lies the usage of the ð¿0-normrnminimization of the objective function, which consists of data misfit and L0-normrnstabilizing function, by an efficient iteratively reweighted least-squares (IRLS)rnalgorithm. As a major contribution, the presented method incorporates three novelrntechnical advancements. At first, the method incorporates an auto-adaptivernregularization technique, which automatically determines a suitable regulariza t ionrnparameter, and a modified error weighting function that helps to improve both thernstability and convergence of the method. The other advantage of the auto-adaptivernregularization technique and error weighting matrix is that they are not whollyrndependent on the known noise level. Because of that, the method can yield reasonablernresults even when the noise level of the data is not properly known.rnThe other major contribution is the use of a new depth weighting function. Thernadvantages of the newly proposed depth-weighting function can be summarized asrnfollows: (I) It properly counteracts the gravity kernel decay, so that the inversion resultsrnrncan provide realistic depth information. (II) It avoids the selection of the depthrnweighting function parameters through trial and error, which was common in therntraditional methods, through the usage of automated parameters selection techniques.rnEspecially, this has a significant advantage when there is no prior depth informa t ionrnthat helps to choose the optimum parameter in using the traditional depth weightingrnfunction.rnFurthermore, to achieve geologically plausible results and also reduce solutionrnambiguity, a physical parameter inequality constraint algorithm has been developed andrnemployed to constrain the obtained density contrast values.rnFinally, the implementation of an effectively combined stopping criterion hasrnbeen used to terminate the iterative inversion procedure, when geologically viablernsolutions are obtained. The proposed combined termination criteria are shown tornoutperform traditional termination criteria, used in most iterative geophysical invers ionrnalgorithms, through the modeling of synthetic and published measured datarnTo test the new method, forward and inverse modeling codes were written usingrnPython programming in a Linux environment. The validity of the overall Python codes,rnand the practicality and efficiency of the presented inversion method were tested byrninverting a number of synthetic data sets from geometrically complex bodies and fieldrndata sets from different geological settings. The results of the inversion confirmed therncapability of the developed inversion method in producing geologically acceptablerncompact and sharp models that define the shape, location, and density distribution ofrnthe causative subsurface bodies. This proves, the reliability and effectiveness of therndeveloped inversion method in practical applications to delineate sharp discontinuit iesrnand blocky features such as distinct layering or formation of localized bodies in thernsubsurface.