Search

Dispersion is an important effect contributing to attenuating ground-penetrating radar (GPR) signals and must be incorporated in wave propagation. The constant-[Formula: see text] model (quality factor [Formula: see text]) is valid in rocks for typical GPR frequency ranges. In time-domain wave modeling, to simulate a constant-[Formula: see text] medium, using multiple Debye poles, the number of poles and their parameters must be previously estimated to honor target values for [Formula: see text] and velocity in the demanded frequency band. For electromagnetic waves, this problem is usually formulated as a highly nonlinear optimization problem presenting several solutions, which is solved using global optimization methods. An efficient linear approach to estimate multiple Debye poles to simulate constant-[Formula: see text] media for GPR data is presented. The resulting inverse problem is easy to solve because the number of poles is relatively small. The approach allows for a novel rescaled solution that depends only on the frequency band and number of poles. Solutions satisfying specific target values for velocity and [Formula: see text] are obtained from the rescaled solution by multiplying this solution by a proportionality factor that depends only on the target values. Choosing the number of poles is done so that velocity [Formula: see text] and [Formula: see text] curves, as a function of frequency [Formula: see text], are approximately invariant in the specified frequency band. Using the criterion that [Formula: see text] and [Formula: see text] curves must present maximum percentage deviations from the respective target values smaller than 5%, we find that, for the 200 MHz antenna, whose useful frequency band is 50–800 MHz, only three poles are sufficient. The theory is validated by modeling dispersion effects using synthetic models. The importance of introducing dispersion effects is illustrated by presenting a trial-and-error modeling approach to approximately reproduce a 200 MHz radargram acquired in fractured and karstified carbonate rocks.

Most of the approaches designed to avoid cycle skipping in full-waveform inversion (FWI) involve calculating a sequence of inversions in a multiscale fashion. We have adopted an alternative strategy, which is inverting a sequence of different misfit functions in the time domain. This is an implicit multiscale approach in the sense that the used misfit functions are sensitive to different wavelengths, but all of the inversion steps use the same modeling algorithm and the same model grid. In the first and third inversion steps, the transmitted (early arrivals) and reflected (late arrivals) components of the wavefield envelopes are respectively fitted. The second step promotes a smooth transition between the first and third steps, by using the envelope of the complete waveform. Because fitting just the envelope of the reflected waves has a minor effect on the misfit function of the whole data set, the phases of the reflected waves are mostly fitted in the fourth step, which is based on the waveform misfit function preserving only the late arrivals. The third and fourth steps are of crucial importance to fit the reflected events. We test the sequential inversion approach with the Marmousi model using data sets with different frequencies, obtaining better estimates of the velocity field than those obtained with the classic FWI. The solutions obtained with classic FWI and sequential inversion approach degrade with a progressively higher peak frequency data set, but the classic FWI solution degrades more rapidly.

We have developed an inversion approach based on particle swarm optimization (PSO) for 2D direct current resistivity data. Once the PSO algorithm is implemented, cost-function changes on constraints and misfit criteria can be easily performed and an approximate invariance of the Lagrange multipliers is obtained if each term of the cost function is properly normalized. More importantly, the interpreter can divide the mesh into partitions and apply a different constraint or combination of constraints in each partition. We explore several restrictions on the log-resistivity variation, either in single or in joint versions, including spatial continuity constraints in the [Formula: see text]-and [Formula: see text]-norms, total variation, sparsity constraints using the discrete cosine transform and Daubechies bases, and minimum moment of inertia constraints to impose concentration of resistive or conductive materials along target axes. In the latter case, the earth surface might be used also as a target axis. In addition, the final state of the PSO algorithm, once a stopping condition has been reached, contains not only the best solution but also a cluster of good suboptimal solutions that can be used for uncertainty analysis. As a result, the interpreter has the flexibility to perform an interpretation process based on a feedback trial-and-error inversion approach, in a similar manner that he/she has when using a friendly forward modeling software. The interpreter can then drive the tentative solutions obtained along the inversion process to incorporate their conceptions about the geologic environment, besides appraising misfit and stability of the solutions. We evaluate synthetic and field data examples. In the first case, we determine how the interpreter can drive the inversion process in a karst environment to image dissolution features. In the second case, we explore a similar situation aiming to optimize borehole locations in fracture zones in crystalline rocks.

One-dimensional forward modeling in direct current (DC) resistivity is actually computationally inexpensive, allowing the use of global optimization methods (GOMs) to solve 1.5D inverse problems with flexibility in constraint incorporation. GOMs can support computational environments for quantitative interpretation in which the comparison of solutions incorporating different constraints is a way to infer characteristics of the actual subsurface resistivity distribution. To this end, the chosen GOM must be robust to changes in the cost function and also be computationally efficient. The performance of the classic versions of the simulated annealing (SA), genetic algorithm (GA), and particle swarm optimization (PSO) methods for solving the 1.5D DC resistivity inverse problem is here compared using synthetic and field data. The main results are as follows: (1) All methods reproduce synthetic models quite well, (2) PSO and GA are comparatively more robust to changes in the cost function than SA, (3) PSO first and GA second present the best computational performances, requiring less forwarding modeling than SA, and (4) GA gives higher performance than PSO and SA with respect to the final attained value of the cost function and its standard deviation. From our experience, to put them into effective operation, the methods can be classified from easy to difficult in the order PSO, GA, and SA as a consequence of robustness to changes in the cost function and of the underlying simplicity of the associated equations. To exemplify a quantitative interpretation using GOMs, we compare solutions with least-absolute and least-squares norms of the discrepancies derived from the lateral continuity constraints of the log-resistivity and layer depth as a manner of detecting faults. GOMs additionally provide the important benefit of furnishing not only the best solution but also a set of suboptimal quasisolutions from which uncertainty analyses can be performed.

Additional information regarding the continuity and resolution of selected seismic reflectors in reverse time migration (RTM) images can be beneficial for seismic interpretation. We have developed and evaluated new imaging conditions for RTM based on the phase coherence between the forward- and backward-propagated wavefields. These imaging conditions make use of the instantaneous phase and envelope of the analytical signals of the source and receiver wavefields, in addition to their real parts. Once the analytical wavefields are available, these imaging conditions can be calculated simultaneously with conventional conditions at little or no extra cost. The availability of these fields at each image point enables several alternative ways to define imaging conditions. We explore, in addition to pure phase crosscorrelation (PC), two approaches of amplitude-weighted PC. Our numerical experiments, imaging synthetic and field data sets, indicate that these new imaging conditions provide additional images that can highlight some weak reflectors by locally improving the resolution of RTM images. In our examples, this happens particularly in the deep portions of the seismic images. In addition, reflection events produced at discontinuities are enhanced as sharp signals, suggesting that the proposed imaging conditions can help to delineate stratigraphic and structural features that are harder to see in conventional images. These properties of the PC imaging conditions make them an interesting tool to provide additional information that can aid seismic interpretation in complex structural settings. © 2018 Society of Exploration Geophysicists., We thank CNPq for the financial support through the Instituto Nacional de Ciências e Tecnologia em Geofísica do Petróleo (INCT-GP). J. C. Costa, W. E. Medeiros, and J. Schleicher also thank CNPq for the research fellowships and associated grants. M. Schimmel thanks the Science without Borders Program for his PVE grant and Misterios project (CGL2013-48601-C2-1-R).

We have developed a new method for the approximate mapping of subsurface resistivity based on a geostatistical approach. We have used the characteristic points (CPs) of a resistivity sounding, which are their inflection and extreme (maximum and minimum) points. The methodology consists of five stages: (1) smoothing the geoelectric soundings to assure stability in the presence of measurement errors, (2) identifying the CPs in the smoothed versions of the geoelectric soundings, (3) obtaining subsurface sparse point resistivity estimates for the subsurface resistivity from the CPs, using empirical relationships between electrode spacing and penetration depth, (4) calculating semivariograms using the subsurface sparse point resistivity estimates, and fitting them with semivariogram models, and, finally, (5) estimating the subsurface resistivity distribution by kriging interpolation of the sparse point estimates. No explicit assumption is made about the true subsurface resistivity and, as result, the method can be applied to models of any dimension, given that the empirical rule for the penetration depth is approximately valid. Computationally, the method is very fast because no numerical modeling (either direct or inverse) is demanded and the most intensive computer operation is just a kriging interpolation. Given its independence of the model dimension and stability in the presence of measurement errors, the method can be readily implemented as a fast automatic method of interpretation. As result, the best scenario to apply the method is the initial treatment of large data sets. The methodology was developed, with synthetic and field data, using direct current (DC) resistivity soundings obtained with the Schlumberger array. The proposed method might be generalized for other DC-resistivity arrays and electromagnetic techniques based on resistivity soundings.

The key aspect limiting resolution in crosswell traveltime tomography is illumination, a well-known result but not well-exemplified. We have revisited resolution in the 2D case using a simple geometric approach based on the angular aperture distribution and the Radon transform properties. We have analytically found that if an isolated interface had dips contained in the angular aperture limits, it could be reconstructed using just one particular projection. By inversion of synthetic data, we found that a slowness field could be approximately reconstructed from a set of projections if the interfaces delimiting the slowness field had dips contained in the available angular apertures. On the one hand, isolated artifacts might be present when the dip is near the illumination limit. On the other hand, in the inverse sense, if an interface is interpretable from a tomogram, there is no guarantee that it corresponds to a true interface. Similarly, if a body is present in the interwell region, it is diffusely imaged, but its interfaces, particularly vertical edges, cannot be resolved and additional artifacts might be present. Again, in the inverse sense, there is no guarantee that an isolated anomaly corresponds to a true anomalous body, because this anomaly could be an artifact. These results are typical of ill-posed inverse problems: an absence of a guarantee of correspondence to the true distribution. The limitations due to illumination may not be solved by the use of constraints. Crosswell tomograms derived with the use of sparsity constraints, using the discrete cosine transform and Daubechies bases, essentially reproduce the same features seen in tomograms obtained with the smoothness constraint. Interpretation must be done taking into consideration a priori information and the particular limitations due to illumination, as we have determined with a real data case.

Improving the attributes of Reverse Time Migration (RTM) images can be beneficial for seismic interpretation. We present and evaluate new image conditions for RTM taking also into account the phase coherence between the forward and backward propagated wavefields. To compute the new image conditions we need to propagate the source and receivers analytical signals; as consequence, we have at each image point source and receivers amplitude and phase fields, besides the real wavefields. These fields provide several alternatives to define imaging conditions. Here we explore one of the possibilities that is to weight the usual cross-correlation term with a power n of the phase cross-correlation. Our numerical experiments, imaging the Marmousi data, show that tunning the phase cross-correlation exponent parameter n in the imaging condition can produce sharp images across layer boundaries, suggesting that phase cross-correlation in RTM helps to enhance signals, reduce backscatering artfacts and highlight subtle stratigraphic and structural features. The properties of the phase cross-correlation imaging conditions can be advantageous for seismic interpretation in complex structural settings.

We illustrate the importance of establishing solution uniqueness through mathematical restrictions reflecting a source attribute. We also illustrate the validity and utility of a guideline derived in an accompanying paper for constructing sound gravity inversion methods for the class of sources presenting either homogeneous or depth‐independent density distributions. The two‐part guideline is (1) to introduce a priori information favoring uniqueness, either by assuming that a nonnull density distribution depending only on x and y is confined to the interior of a horizontal slab with known position or by limiting the class of possible solutions to homogeneous, simply connected polygons (or polyhedra) with known density, displaying no fancy shapes and no curling apophyses at their borders, and (2) to introduce information favoring solution stability by estimating only the features of the source which may be resolved by the data. Following the guideline, we apply different methods to gravity data using interpretation models consisting of a grid of cells on the x‐y and x‐z planes. In both cases the estimates are very close to the true synthetic source. The data produced by the distribution varying with x and z are also inverted using the method, which minimizes the norm of the first‐order derivative of the density. This constraint does not reflect a true source attribute but is strong enough to stabilize the solution and to guarantee its uniqueness. Because of the strong bias imposed to the solution, the estimated distribution, although unique and stable, is far from the true source, concentrating most of the anomalous mass at the surface. Finally, we present an alternative method which redistributes the estimated anomalous mass downward. To be effective, this technique requires prior knowledge about the source depth to the top. In addition, the source should not be too small and deep. Although being able to produce good results, this alternative method requires a great dose of the interpreter's art.

To obtain a unique and stable solution to the gravity inverse problem, a priori information reflecting geological attributes of the gravity source must be used. Mathematical conditions to obtain stable solutions are established in Tikhonov's regularization method, where the a priori information is introduced via a stabilizing functional, which may be suitably designed to incorporate some relevant geological information. However, there is no unifying approach establishing general uniqueness conditions for a gravity inverse problem. Rather, there are many theorems, usually establishing just abstract mathematical conditions and making it difficult to devise the type of geological information needed to guarantee a unique solution. In Part I of these companion papers, we show that translating the mathematical uniqueness conditions into geological constraints is an important step not only in establishing the type of geological setting where a particular method may be applied but also in designing new gravity inversion methods. As an example, we analyze three uniqueness theorems in gravimetry restricted to the class of homogeneous bodies with known density and show that the uniqueness conditions established by them are more probably met if the solution is constrained to be a compact body without curled protrusions at their borders. These conditions, together with stabilizing conditions (assuming a simple shape for the source), form a guideline to construct sound gravity inversion methods. A historical review of the gravity interpretation methods shows that several methods implicitly follow this guideline. In Part II we use synthetic examples to illustrate the theoretical results derived in Part I. We also illustrate that the presented guideline is not the only way to design sound inversion methods for the class of homogeneous bodies. We present an alternative approach which produces good results but whose design requires a good dose of the interpreter's art.

Despite being widely used, Euler deconvolution has received little attention in its theoretical aspects. We obtain analytical estimators [Formula: see text] and [Formula: see text] for the horizontal [Formula: see text] and the vertical [Formula: see text] source positions, which are the parameters to be estimated when the structural index is assigned a tentative value. By analyzing the estimators’ properties, we show two things. (1) Scattering of the alternative solutions associated with different data window positions has two components: data noise and the use of a wrong tentative value for the structural index, (2) The differences [Formula: see text] and [Formula: see text], as a function of the data window position, exhibit antisymmetric and symmetric behaviors, respectively, about [Formula: see text] when the magnetization and the geomagnetic field inclinations are 90°, 0°, or 45°. For intermediate inclinations there is a moderate departure from the symmetry properties. This analysis shows that the criterion for determining the structural index as the tentative value producing the smallest solution scattering is theoretically sound but occasionally fails in practice because data noise also contributes to solution scattering. In addition, the lower sensitivity of the averages of estimates of the horizontal position (as compared with the vertical position) occurs because [Formula: see text] displays an antisymmetric behavior with respect to the true horizontal position while the estimates of the vertical position display a symmetric behavior with respect to the same point.

To produce a unique and stable solution in potential‐field interpretation, an inversion method must introduce particular constraints. These constraints will inevitably restrict the type of geological setting where the method may be applied. We present a nonmathematical overview of most stabilizing constraints used in inversion methods. Our purpose is to demonstrate that the inversion results are valuable only if the mathematical stabilizing constraints are translated from the geological setting. We identify five basic types of constraints: (1) lower and upper bounds of parameter estimates; (2) proximity of a parameter estimate to a specified value; (3) proximity between pairs of parameter estimates; (4) concentration of the anomalous source about a geometrical element such as an axis; and (5) source compactness. In practice, if used in isolation, constraints (1), (2), (4), and (5) will not produce geologically meaningful results, regardless of the geological setting of the interpretation area. Constraint (3) may produce geologically meaningful results if the anomalous source has a spatially smooth attribute such as the physical property. We illustrate that constraints 1–4, if used in isolation, cannot delineate the geometry of a simulated sill intruded into a sedimentary basin. The basic constraints may (and should) be combined in inversion to produce geologically meaningful results. We present two examples of effective constraint combination: (1) proximity to a specific value and mass concentration about an axis (used to delineate the thickness variation of a sill intruded in a sedimentary basin) and (2) inequality, proximity of a parameter estimate to a specified value, and proximity between pairs of parameter estimates (used to map a discontinuous basement relief). Usually, the stabilizing constraints are too restrictive to hold at all points of a given geological environment. In this case, we use different constraints in different sub‐areas. Each constraint is based on its compatibility with the actual geology of the subarea.

We present a constraint which incorporates in the geophysical inverse problem a priori information about the source convexity. Two kinds of convexity are considered: directional convexity, defined as the attribute of a body being intersected at most at two points by a straight line with a fixed spatial orientation, and global convexity, defined as the attribute of a body being intersected at most at two points by any arbitrarily oriented straight line. The interpretation model consists of several juxtaposed 2-D prisms whose positions, widths, and physical property are established by the interpreter. The depth to the top and the thickness of each prism are the parameters to be determined. The directional convexity is incorporated by minimizing a functional which expresses the number of times that several lines, parallel to a given direction and spaced at regular intervals, intersect the interpretation model. On the other hand, the global convexity is introduced in an algorithmic way by constraining the depths to the top and to the bottom of each prism to be, respectively, smaller or greater than the average of the depths to the top or to the bottom of the adjacent prisms. The use of the presented convexity contraint is demonstrated using synthetic gravity data produced by an isolated source. We found that convexity constraints (either directional along the horizontal axis or global) alone are insufficient to produce strongly stable solutions, although they reduce substantially the solution instability. Despite the “slight” instability, the directional convexity constraint leads to a reasonable delineation of the source shape and produces a better resolution of its base as compared with the inverse method imposing smoothness constraints on the top and bottom surfaces of the anomalous source. The potential use of the convexity constraint is the possibility of combining it with other constraints to produce geologically meaningful solutions. As an example, the global convexity was combined with the minimum moment of inertia of the anomalous masses about a vertical axis passing through its center of mass. This combination allowed introducing a priori information about a diapir shape. This is important because a maximum spread at the diapir middle portion indicates that it ascended like a viscous bubble in an isotropic medium. On the other hand, a spread located at its top indicates that the host medium is mechanically anisotropic, exhibiting horizontal planes of weakness where the mobile material was forced into. So, this combination of constraints alows introducing a priori information about the diapir emplacement. The combination of global convexity and minimum moment of inertia may also be applied to isolated intrusions, as illustrated with the real gravity anomaly produced by the granitic body of Castelsarrasin, in the Aquitaine Basin, France. The estimated source has its top located between 0.45 and 0.57 km from the surface, a base 6.4 km deep, and center of mass at a depth of 3.5 km, which agree reasonably with a priori geological and geophysical information about the body.

We present a new stable gravity inversion method applied to the mapping of an interface separating two homogeneous media. In contrast with previous similar methods, it does not impose an overall smoothness on the estimated interface to stabilize the solution. The density contrast between the media is assumed to be known. The interpretation model for the upper medium consists of rectangular juxtaposed prisms whose thicknesses represent the depths to the interface and are the parameters to be estimated. The true interface is assumed to be flat everywhere except at faults. To incorporate this attribute into the estimated relief, we developed an iterative process in which three kinds of constraints are imposed on parameters: (1) proximity between values of adjacent parameters, (2) lower and upper bounds to parameters, and (3) proximity between the values of parameters and fixed numerical values. Starting with an initial solution which presents an overall smooth relief, the method enhances initially estimated geometric features of the interface; that is, flat areas will tend to become flatter and steep areas will tend to become steeper. This is accomplished by weighting the constraints, which requires proximity between adjacent parameters. The weights are initialized with values equal to unity and are updated automatically to enhance any discrepancy between adjacent depths that have been detected at the initial solution. Constraints 2 and 3 are used both to compensate for the decrease in solution stability caused by the introduction of small weights and to reinforce flatness at the basin bottom. Constraint 2 imposes that any depth be nonnegative and smaller than an a priori known maximum depth value, whereas constraint 3 imposes that all depths be closest to a value greater than the maximum depth. The trade‐off between these conflicting constraints is attained with a final relief presenting flat bottom and steep borders. The method was tested with a synthetic gravity anomaly produced by a simulated sedimentary cratonic extensional basin whose basement consists of steep edges and a flat bottom. The results showed an improvement in the resolution of the relief, leading to a reliable mapping both of the sharp discontinuities at the borders and of the lateral extent of the base of the basin. Additionally, the method produced excellent estimates for the average dip angles of the basin edges (presumably controlled by normal faults), indicating, in this way, its potential in interpreting data produced by this kind of basin. The method was applied to the Bouguer anomaly from the northern portion of Steptoe Valley, Nevada, delineating an isolated basin with a wider, flat base and relatively straight borders as compared with the estimate imposing overall smoothness on the relief.

We present a new, stable method for interpreting the basement relief of a sedimentary basin which delineates sharp discontinuities in the basement relief and incorporates any law known a priori for the spatial variation of the density contrast. The subsurface region containing the basin is discretized into a grid of juxtaposed elementary prisms whose density contrasts are the parameters to be estimated. Any vertical line must intersect the basement relief only once, and the mass deficiency must be concentrated near the earth’s surface, subject to the observed gravity anomaly being fitted within the experimental errors. In addition, upper and lower bounds on the density contrast of each prism are introduced a priori (one of the bounds being zero), and the method assigns to each elementary prism a density contrast which is close to either bound. The basement relief is therefore delineated by the contact between the prisms with null and nonnull estimated density contrasts, the latter occupying the upper part of the discretized region. The method is stabilized by introducing constraints favoring solutions having the attributes (shared by most sedimentary basins) of being an isolated compact source with lateral borders dipping either vertically or toward the basin center and having horizontal dimensions much greater than its largest vertical dimension. Arbitrary laws of spatial variations of the density contrast, if known a priori, may be incorporated into the problem by assigning suitable values to the nonnull bound of each prism. The proposed method differs from previous stable methods by using no smoothness constraint on the interface to be estimated. As a result, it may be applied not only to intracratonic sag basins where the basement relief is essentially smooth but also to rift basins whose basements present discontinuities caused by faults. The method’s utility in mapping such basements was demonstrated in tests using synthetic data produced by simulated rift basins. The method mapped with good precision a sequence of step faults which are close to each other and present small vertical slips, a feature particularly difficult to detect from gravity data only. The method was also able to map isolated discontinuities with large vertical throw. The method was applied to the gravity data from Reco⁁ncavo basin, Brazil. The results showed close agreement with known geological structures of the basin. It also demonstrated the method’s ability to map a sequence of alternating terraces and structural lows that could not be detected just by inspecting the gravity anomaly. To demostrate the method’s flexibility in incorporating any a priori knowledge about the density contrast variation, it was applied to the Bouguer anomaly over the San Jacinto Graben, California. Two different exponential laws for the decrease of density contrast with depth were used, leading to estimated maximum depths between 2.2 and 2.4 km.

Euler deconvolution has been widely used in automatic aeromagnetic interpretations because it requires no prior knowledge of the source magnetization direction and assumes no particular interpretation model, provided the structural index defining the anomaly falloff rate related to the nature of the magnetic source, is determined in advance. Estimating the correct structural index and electing optimum criteria for selecting candidate solutions are two fundamental requirements for a successful application of this method. We present a new criterion for determining the structural index. This criterion is based on the correlation between the total‐field anomaly and the estimates of an unknown base level. These estimates are obtained for each position of a moving data window along the observed profile and for several tentative values for the structural index. The tentative value for the structural index producing the smallest correlation is the best estimate of the correct structural index. We also propose a new criterion to select the best solutions from a set of previously computed candidate solutions, each one associated with a particular position of the moving data window. A current criterion is to select only those candidates producing a standard deviation for the vertical position of the source smaller than a threshold value. We propose that in addition to this criterion, only those candidates producing the best fit to the known quantities (combinations of anomaly and its gradients) be selected. The proposed modifications to Euler deconvolution can be implemented easily in an automated algorithm for locating the source position. The above results are grounded on a theoretical uniqueness and stability analysis, also presented in this paper, for the joint estimation of the source position, the base level, and the structural index in Euler deconvolution. This analysis also reveals that the vertical position and the structural index of the source cannot be estimated simultaneously because they are linearly dependent; the horizontal position and the structural index, on the other hand, are linearly independent. For a known structural index, estimates of both horizontal and vertical positions are unique and stable regardless of the value of the structural index. If this value is not too small, estimates of the base level for the total field are stable as well. The proposed modifications to Euler deconvolution were tested both on synthetic and real magnetic data. In the case of synthetic data, the proposed criterion always detected the correct structural index and good estimates of the source position were obtained, suggesting the present theoretical analysis may lead to a substantial enhancement in practical applications of Euler deconvolution. In the case of practical data (vertical component anomaly over an iron deposit in the Kursk district, Russia), the estimated structural index (corresponding to a vertical prism) was in accordance with the known geology of the deposit, and the estimates of the depth and horizontal position of the source compared favorably with results reported in the literature.

We present a gravity interpretation method for estimating the relief of an arbitrary interface separating two homogeneous media. The upper medium is discretized into rectangular, juxtaposed prisms whose thicknesses represent the depths to the interface and are the parameters to be estimated from the gravity anomaly. The density contrast of each prism is assumed to be constant and known. To stabilize the inversion, we introduce two kinds of constraints on the depths. The first one requires proximity between the observed and computed depths at isolated points such as those obtained from boreholes (absolute equality constraint). The second one requires that groups of depths approximately follow an established linear relationship among the depths (relative equality constraint). Both kinds of constraints are imposed in the least‐squares sense. We illustrate the method performance by applying it to a synthetic anomaly produced by a simulated basement relief consisting of four narrow and adjacent structural lows. Only two structural lows produced isolated gravity lows. Nonetheless, the whole basement topography was successfully reconstructed with an average error of 4% of the maximum relief amplitude. In this example, the relative constraints established that the thicknesses of adjacent prisms should be as close as possible to each other (overall relief smoothness). As absolute constraints we used point information about the basement depth at five points. In addition to the overall relief smoothness, other relevant geologic information, such as localized relief smoothness, occurring at structural terraces, for example, can be incorporated by assigning different weights to the relative equality constraints. The method was applied to the gravity anomaly of Recôncavo Basin, Brazil, leading to a sharper definition of the basement features relative to previous gravity interpretations of this area.

Current inversion methods estimating a physical property on a grid use two kinds of constraints: (1) absolute constraints requiring that the solution vector be close to a prespecified vector as, for example, proximity to the null vector in the ridge regression, and (2) relative constraints requiring that the elements of the solution vector be close to each other as, for example, smoothness imposed by methods requiring that spatial derivatives of the physical property be continuous. Using just absolute constraints may produce a biased solution conflicting with some true geological attributes of the source. On the other hand, using just relative constraints may be insufficient to stabilize the inversion. We present a stable inversion method incorporating absolute constraints only at points where the physical property is known as, for example, outcroppings and boreholes. Elsewhere, relative constraints are introduced according to a linear relationship for the spatial distribution of the physical property. The proposed method is analogous to an interpolation method applied to the physical property distribution: the interpolating function must satisfy some property, like continuity, and pass through the data. The proposed method is versatile in incorporating not only factual, but also indirect information such as continuity, symmetry, and trends. In addition, it may be applied both to linear and nonlinear problems. The relevant results obtained from applying the method to synthetic potential‐field data were: (1) ability to produce reliable mappings of the basement relief of a sedimentary basin (from gravity data) and of an areal distribution of magnetization, using only smoothness constraints and restricted surface and borehole information; (2) operational simplicity because of the solution insensitivity to the damping parameter. The method is also applied to the Bouguer anomaly over San Jacinto Graben, California, producing results consistent with previous interpretations.

Magnetic interpretations are usually carried out either by assuming induced magnetization and estimating the model geometry, or by presuming a known source spatial orientation to estimate the total magnetization. We present a 3-D magnetic interpretation method that estimates simultaneously the total magnetization direction and the spatial orientation of the source. It is based on the approximation of the anomaly by the series derived from expanding the magnetic potential into multipoles and retaining source moments up to second order. The moments and linear combinations of moments appearing in the series are then inverted from the magnetic anomaly. The total magnetization is assumed constant in direction but not in magnitude. It is also presumed implicitly that the anomalous distribution of magnetization intensity has nonzero values in a finite‐volume region, is far from the observation points, and presents three othogonal planes of symmetry intersecting at the center of the dipole moment. The method is essentially linear and requires no a priori explicit assumption of a fixed geometry for the sources. The method is particularly suited to interpret compact, isolated or disjoint, but spatially correlated sources. If the source satisfies all assumptions presumed by the method, it is possible to obtain accurate, stable estimates of the total dipole moment vector, the position of the center of dipole moment, and the directions of all three principal axes of symmetry. If the source is not far from the observation plane and/or if the total magnetization direction is not constant, it is still possible to obtain accurate and stable estimates of the direction of the mean total magnetization and the projection, on the observation plane, of the center of dipole moment. The method is applied to magnetic data from the Gulf of Guinea Seamount. The estimated magnetic palaeopole is at 50°48′S and 74°54′E which is in good agreement with estimates published by other authors.

We present a 3-D gravity interpretation method based on the inversion of source moments employing the series derived from the gravity anomaly expansion in multipoles and retaining moments up to second order only. It presents the advantages of being essentially linear and allowing a straightforward computer implementation. The method requires neither an explicit assumption about the source geometry nor a priori knowledge about the density contrast distribution, which may even be nonuniform. The method assumes implicitly: (1) that the source fits in the interior of an imaginary sphere whose center coincides with the source center of mass and whose radius is smaller than the depth to the source center of mass relative to the measuring plane, and (2) that the spatial density distribution presents three ortogonal planes of symmetry intersecting at the source center of mass. The first assumption can be met by upward continuing the observed anomaly. When both assumptions are met, the method produces reliable and stable estimates of the total anomalous mass, the coordinates of the center of mass, the three principal axes directions of the anomalous body, and the relative importance among the axes. The method is particularly suited for interpreting compact, isolated or disjoint, but spatially correlated sources. The method is sensitive to an incomplete removal of the regional field (presumably superimposed on the residual anomaly of interest) because the unremoved part of the regional field may strongly affect the estimates of the residual source moments.

Euler deconvolution has been widely used in automatic aeromagnetic interpretations because it requires no knowledge of the source magnetization direction and assumes no particular interpretation model. This rapid inversion method produces stable estimates for the horizontal and vertical positions and respectively) of point sources. Another main advantage of this method is that it does not demand a detailed geologic knowledge of the area. However the two crucial weakenesses of this method are: i) the need of knowing the nature of the magnetic point or line source, that is, the knowledge of the structural index and ii) the presence of an undesirable “spray” of solutions which is obtained when the method is applied in a redundant way by moving the discrete inverse operator around the anomalous area. We present two new criteria for solving these difficulties. The first criterion is aimed at determining It is based on the correlation between the total-field anomaly and the estimates of a baselevel (b) computed at the same locations of the observations and using different This correlation attains a minimum when the correct is assumed. The second criterion selects, from the “spray” of solutions, only the best solutions. We propose that, besides selecting the solutions producing a standard deviation for smaller than a threshold value, only those solutions producing the “best fit” to the known quantities be selected. This additional criterion substantially reduces the “spray” of solutions.

We present a gravity interpretation method for estimating the relief of an arbitrary interface separating two homogeneous media. The upper medium is discretized into rectangular, juxtaposed prisms whose thicknesses represent the depths to the interface and are the parameters to be estimated from the gravity anomaly. To stabilize the inversion, we introduce two kinds of constraints on the depths. The first one requires proximity between the observed and computed depths at isolated points such as those obtained from boreholes (absolute equality constraint). The second one requires that groups of depths approximately follow an established linear relationship among them (relative equality constraint). Both kinds of constraints are imposed in the least-squares sense. This stabilization procedure does not bias the solution either toward an initial approximation, or toward the null solution. The method is applied to synthetic data, producing a good reconstruction of a complex basement consisting of several structural high and lows. It is also applied to real data of Reconcavo Basin, Brazil, leading to a sharper definition of the basement structures as compared with previous gravity interpretations of the area.

We present theoretical uniqueness, stability, and parameter bias analyses for the simultaneous estimation of: the structural index the horizontal and vertical positions and respectively) of a point magnetic source, and a baselevel b, using Euler’s differential equation. We found that and are independent of the other parameters and can be estimated from the knowledge of the total-field anomaly and its gradients only. On the other hand, and cannot be simultaneously estimated because they are linearly dependent. As a result, to estimate must be either known a priori or independently determined (and vice versa). For a known estimates of and are unique and stable regardless of the value of However, if the vertical (or horizontal) gradient of the observed total field assumes negligible values, the estimation of (or becomes unstable. In contrast, as far as is not too close to zero, estimates of are stable. For equal to zero, estimates of are undetermined.

The interpretation of magnetic anomalies is not a trivial task not only because of the difficulties in allowing for a possible existence of remanent magnetization, but also due to the nonlinear dependence of the magnetic anomaly upon the source geometry. For the sake of simplicity, the magnetization is usually assumed to be induced only while trying to estimate the interpretation model geometry. Alternatively, the spatial orientation of the source may be assumed known in order to estimate the total magnetization. We present a method to simultaneously estimate the total magnetization direction and the spatial orientation of the source. The method is based on the approximation of the anomalous field by the series derived from expansion of the magnetic potential in multipoles, keeping moments up to 2nd order. Because the parameters to be inverted are the source moments, the method has the advantages of being linear and appropriate to interpret anomalies produced by a broad class of sources. This class comprises the finite-volume bodies, homogeneously magnetized in direction, but not necessarily in modulus, and whose spatial distribution of magnetization modulus presents three orthogonal symmetry planes.

Gravity inversion is an ill-posed problem because it presents infinite solutions. One efective way to transform it into an well-posed problem is to constrain the solutions with a priori information (about the anomalous density distribution) introduced by a stabilizing functional. In this paper we present a stabilizing functional which incorporates information about symmetry and either overall or directional smoothness in the spatial density distribution. In the case of a linear inverse problem where the parameters to be determined are the densities of several retangular blocks, the proposed functional is the sum of the squared differences between previously established pairs of parameters. The estimator associated with this functional may contribute to stabilize the inversion process if it incorporates information about the density variation in a direction in which the gravity observations have poor resolution particularly in the vertical direction. Let us assume that by inspecting a given gravity profile, we can limit a region R containing the cross-section l? of the anomalous source. To estimate, from the gravity profile, the distribution p(z, z) of density contrasts, we discretize R into M blocks (Figure 1) and estimate them by solving the equation