The prism method for terrain corrections using digital computers

Summary In the prism method of making terrain corrections, the topography is approximated by a model consisting of right rectangular prisms. The vertical component of the gravitational attraction of each prism is calculated and the sum of these components gives the terrain correction.The prism method as programmed has no computational limitations. It can be used on all sizes of computers; it can be applied to a large area with any fine grid interval; it can be processed in a single run and yet provides complete flexibility for both research and routine computations. This has been achieved by breaking up larger areas into regions which fit into the computer memory. The contributions of these regions are automatically summed up for each station. While processing each region, various controls may be used at each station to exclude the contribution of a distant part of the area, to use approximate expressions farther from the station, to print out details around the station. There is also provision to refine the model by using smaller prisms around each computation point. Thus full use of elevation control can be made to calculate the terrain correction, the accuracy of which depends only on the quality of the input data.The prism method has been used to calculate terrain corrections for 130 stations in the New Quebec crater area. For five of these stations terrain corrections were also calculated by usingHammer's template. The two independent sets of values differ by less than four per cent.     

Improving the computation of the gravitational terrain effect close to ground stations in the GTE software

The precise computation of the vertical gravitational attraction of the topographic masses (terrain correction) is still being studied both for geodetic and geophysical applications. In fact, it is essential in high precision geoid estimation by means of the well-known remove-compute-restore technique, which is used to isolate the gravitational effects of anomalous masses in exploration geophysics. The terrain correction can be evaluated exploiting a Digital Terrain Model (DTM) in different ways, such as classical numerical integration, prisms, tesseroids, polyhedrons, and/or Fast Fourier Transform techniques. The increasing resolution of recently developed DTMs, the increasing number of observation points, and the increasing accuracy of gravity data represent, nowadays, major challenges for the terrain correction computation. Classical point mass approximation and prism based-algorithms are indeed too slow, while Fourier-based algorithms are usually too much approximate when compared to the required accuracy. In this work, we improve the Gravity Terrain Effects (GTE) algorithm, the innovative tool that exploits a combined prism-Fast Fourier Transform approach especially developed for airborne gravimetry, to compute the terrain correction on the surface of the DTM (i.e. corresponding to the ground stations and/or its vicinity). This required development of a proper adjustment of the algorithms implemented within the GTE software and also to define and implement a procedure to overcome the problems of the computation of the gravitational effects due to the actual slope of the terrain close to the stations. The latter problem is thoroughly discussed and solved by testing different solutions like concentric cylindrical rings, triangulated polyhedrons, or ultra-high resolution squared prisms. Finally, numerical tests to prove the temporal efficiency and the computational performances of the improved GTE software to compute terrain correction for ground stations are also presented.     

GRAPHICAL EVALUATION OF MAGNETIC AND GRAVITY ATTRACTION OF THREE-DIMENSIONAL BODIES*

A nomogram is presented which enables evaluation of the components of magnetic attraction of a homogeneous finite rectangular prism, and of gravitational attraction due to a uniform rectangular lamina. In practice any three-dimensional body could be approximated by a number of right rectangular prisms of varying dimensions governed by the shape of the body. The magnetic attraction of the whole body is then obtained by numerical summation of the effects of the constituent prisms. For evaluating the gravitational effect, the cross-section of the body corresponding to each elevation contour is approximated by a number of rectangular laminae (or by a stepping polygon) the attraction of which can be determined with the aid of the same nomogram. The total gravitational attraction of the body is obtained by a process of graphical integration along the vertical axis.     

SOME COMMENTS ON THE CALCULATION OF THE GRAVITATIONAL AND MAGNETIC ATTRACTION OF A HOMOGENEOUS RECTANGULAR PRISM *

The use of arctangents rather than arcsines in the expression for the gravitational attraction of a homogeneous rectangular prism reduces computational difficulties. Once a subroutine is available to compute one component of attraction in a Cartesian coordinate system, the other components may be obtained by cyclic permutation of the field point and body coordinate parameters. This technique also readily provides derivatives of the gravitational attraction and hence forms a compact method for the calculation of a magnetic anomaly due to a homogeneous rectangular magnetic prism.     

Analytic Expressions for the Gravity Gradient Tensor of 3D Prisms with Depth-Dependent Density

Variable-density sources have been paid more attention in gravity modeling. We conduct the computation of gravity gradient tensor of given mass sources with variable density in this paper. 3D rectangular prisms, as simple building blocks, can be used to approximate well 3D irregular-shaped sources. A polynomial function of depth can represent flexibly the complicated density variations in each prism. Hence, we derive the analytic expressions in closed form for computing all components of the gravity gradient tensor due to a 3D right rectangular prism with an arbitrary-order polynomial density function of depth. The singularity of the expressions is analyzed. The singular points distribute at the corners of the prism or on some of the lines through the edges of the prism in the lower semi-space containing the prism. The expressions are validated, and their numerical stability is also evaluated through numerical tests. The numerical examples with variable-density prism and basin models show that the expressions within their range of numerical stability are superior in computational accuracy and efficiency to the common solution that sums up the effects of a collection of uniform subprisms, and provide an effective method for computing gravity gradient tensor of 3D irregular-shaped sources with complicated density variation. In addition, the tensor computed with variable density is different in magnitude from that with constant density. It demonstrates the importance of the gravity gradient tensor modeling with variable density.     

Nomograms for determining the gravity attraction of a right rectangular prism

Summary A rectangular prism is one of the versatile models which is often used to approximate many geological features and also for calculation of topographic and isostatic effects. The analytical formula for the gravity attraction of a right rectangular prism contains 24 lengthy expressions which though compatible for use on an electronic computer, are rather prohibitive for hand computations.Some need has, therefore, been felt for designing a simple graphical device to enable a fairly quick estimate of the gravity effect of a rectangular block model, especially for small-scale jobs involving only a few computation points. In so far as is known to the author, no such graphical device has been made available to the geophysicists.The set of nomograms presented in this communication should fulfill this long felt need. With this aid now it is possible to compute in only a few minutes the gravity attraction of a finite rectangular prism at a point lying outside or on the prism surface.

---

Zusammenfassung Das Quader ist eines der Modelle, das häufig als Annäherung vieler geologischer Körper und bei der Berechnung topographischer und isostatischer Effekte angewandt wird. Die analytische Formel für die Schweresanziehung eines aufrechten Quaders enthält 24 lange Formeln, diezwar für elektronische Rechenautomat aber nicht für Berechnungen von Hand geeignet sind.Es bestand deshalb ein Bedarf für ein einfaches graphisches Mittel, welch eine rasche Abschätzung des Schwereseffektes eines Block-models ermöglichen kann, insbesondere für kleinere Aufgaben mit nur wenigen Berechnungspunkten. Dem Verfasser sind noch keine solche Diagramme bis jetzt bekannt.Die Nomogramme dieser Abhandlung sollten diesen Bedarf erfüllen; an Hand dieses Hilfsmittels sollte es möglich sein in wenigen Minuten die Schweresanziehung für jeden Punkt ausserhalb oder auf der Oberfläche eines begrenzten Quaders zu berechnen.

---

---

Paper presented at the 32nd E.A.E.G. Meeting in Edinburgh (May 1970).     

分段光滑曲线边界波动方程数值模拟研究

矩形网格有限差分法在地震波传播数值模拟方面具有计算速度快的显著优势,但该方法在处理复杂边界问题上存在着效率低的严重缺陷.本文针对分段光滑曲线边界定义了尖点处的一种正则导数,给出了矩形网格情形分段光滑曲线网格边界点法向导数的一种插值计算方法.采用矩形网格有限差分法对复杂边界地球介质模型进行地震波场数值模拟,并采用波场系列快照技术揭示地震波在起伏地表和复杂介质中的传播规律.模拟结果表明：法向导数插值计算方法为矩形网格有限差分法处理复杂边界提供了有效途径,采用波场系列快照技术可以清晰地展现地震波在反射界面的反射和透射规律、在尖点的绕射规律以及在自由表面的直达波和多次反射规律.     

The gravitational potential and its derivatives of a right rectangular prism with depth-dependent density following an n-th degree polynomial

The direct gravity problem and its solution belong to the basis of the gravimetry. The solutions of this problem are well known for wide class of the source bodies with the constant density contrast. The non-uniform density approximation leads to the relatively complicated mathematical formalism. The analytical solutions for this type of sources are rare and currently these bodies are very useful in the gravimetrical modeling. The solution for the vertical component of the gravitational attraction vector for the 3D right rectangular prism is known in the geophysical literature for the density variations described by the 3-rd degree polynomial. We generalized this solution for an n-th degree, not only for the vertical component, but for the horizontal components, the second-order derivatives and the potential as well. The 2D modifications of all given formulae are presented, too. The presented general solutions, which involve a hypergeometric functions, can be used as they are, or as an auxiliary tool to derive desired solution for the given degree of the density polynomial as a sum of the elementary functions. The pros-and-cons of these approaches (the complexity of the programming codes, runtimes) are discussed, too.     

TOPOGRAPHIC REDUCTION OF GRAVITY MEASUREMENTS BY NUMERICAL INTEGRATION OF BOUNDARY INTEGRALS*

A numerical method for calculating the topographic reduction of gravity measurements is developed which follows the approximation of the topography by a single valued function. The method involves the conversion of the volume integral for the gravity effect into a two-dimensional definite integral. The definite integral is partly solved by explicit, and partly by numerical, integration, using the Gauss-Legendre quadrature formula. This method is well suited to calculating the topographic reduction of 50 to approximately 1000 m from the station–especially for microgravimetric surveys in areas of steeply sloping terrain. To test the method in practice it was applied in an area of rough relief in Keban (East Turkey).     

Effects of terrain smoothing on topographic shielding correction factors for cosmogenic nuclide‐derived estimates of basin‐averaged denudation rates

Estimation of spatially averaged denudation rates from cosmogenic nuclide concentrations in sediments depends on the surface production rates, the scaling methods of cosmic ray intensities, and the correction algorithms for skyline, snow and vegetation shielding used to calculate terrestrial cosmogenic nuclide production. While the calculation of surface nuclide production and application of latitude, altitude and palaeointensity scaling algorithms are subjects of active research, the importance of additional correction for shielding by topographic obstructions, snow and vegetation is the subject of ongoing debate. The derivation of an additional correction factor for skyline shielding for large areas is still problematic. One important issue that has yet to be addressed is the effect of the accuracy and resolution of terrain representation by a digital elevation model (DEM) on topographic shielding correction factors. Topographic metrics scale with the resolution of the elevation data, and terrain smoothing has a potentially large effect on the correction of terrestrial cosmogenic nuclide production rates for skyline shielding. For rough, high‐relief landscapes, the effect of terrain smoothing can easily exceed analytical errors, and should be taken into account. Here we demonstrate the effect of terrain smoothing on topographic shielding correction factors for various topographic settings, and introduce an empirical model for the estimation of topographic shielding factors based on landscape metrics. Copyright © 2008 John Wiley and Sons, Ltd.     

Gravimetric survey terrain correction using linear analytical approximation

Various methods for computing the terrain correction in a high‐precision gravity survey are currently available. The present paper suggests a new method that uses linear analytical terrain approximations. In this method, digital terrain models for the near‐station topographic masses are obtained by vectorizing scan images of large‐scaled topographic maps, and the terrain correction computation is carried out using a Fourier series approximation of discrete height values. Distant topography data are represented with the help of digital GTOPO30 and Shuttle Radar Topography Mission cartographic information. We formulate linear analytical approximations of terrain corrections for the whole region using harmonic functions as the basis of our computational algorithm. Stochastic modelling allows effective assessment of the accuracy of terrain correction computation. The Perm Krai case study has shown that our method makes full use of all the terrain data available from topographic maps and digital terrain models and delivers a digital terrain correction computed to a priori precision. Our computer methodology can be successfully applied for the terrain correction computation in different survey areas.     

Three-Dimensional Gravity Modeling In All Space

We review available analytical algorithms for the gravity effect and gravity gradients especially the vertical gravity gradient due to a right rectangular prism, a right polygonal prism, and a polyhedron. The emphasis is placed on an investigation of validity, consistency, and especially singularities of different algorithms, which have been traditionally proposed for calculation of the gravity effect on ground (or outside anomalous bodies), when they are applied to all points in space. The rounding error due to the computer floating point precision is estimated. The gravity effect and vertical gradient of gravity in three dimensions caused by a cubic model are calculated by different types of algorithms. The reliability of algorithms for the calculation of gravity of a right polygonal prism and a polyhedron is further verified by using a regular polygonal prism approximating a vertical cylinder and a regular polyhedron approximating a sphere, respectively. By highlighting Haáz-Jung-Plouff and Okabe-Steiner-Zilahi-Sebess' formulae for a right rectangular prism, Plouff's algorithm for a right polygonal prism, and Gouml;tze and Lahmeyer's algorithm for a polyhedron and removing their singularities, we demonstrate that these formulae and algorithms can be used to model the gravity anomaly and its vertical gradient at all possible computation positions.     

RAPID COMPUTATION OF THE GRAVITATION ATTRACTION OF TOPOGRAPHY ON A SPHERICAL EARTH*

A brief review of the existing methods of gravity reduction is given and a new method suitable for use on high speed digital computers is described. The method is based on the formula for the gravitational attraction of a frustum of a cone. The topographic contours are represented by polygons and the x and y coordinates of corners of the polygons constitute the input to the computer. The vertical component of the gravitational attraction is calculated by evaluating the cone formula for a number of vertical sections of the topography. Each vertical section is simplified by adopting a procedure of grouping and averaging for the distant points of the section. The effect of the earth's sphericity is taken into account by lowering the distant points of the sections by amounts determined by the curvature. The computations include the area close to the point at which the attraction is required and may be limited to an area defined by a circle centered at this point. The method is therefore compatible with the conventional zone chart methods. As an illustration of the method the gravitational attraction of Caryn Seamount in the Atlantic Ocean is computed. The total Bouguer correction and the Terrain correction are also computed for an area in northwestern South America and comparisons are made with hand computations by a zone chart method. As an example, for work at sea, the Bouguer corrections for an area near the Island of Mauritius in the Indian Ocean are computed and the effects of sphericity and three-dimensionality are calculated. The gravitational attraction of two-dimensional bodies can be computed in a very similar manner. The attraction of the Puerto Rico Trench model is computed and the results are compared with other methods. The effects of sphericity and assumptions involved in extending the models to infinity are discussed.     

Regional gravity field modeling based on rectangular harmonic analysis

Regional gravity field modeling with high-precision and high-resolution is one of the most important scientific objectives in geodesy,and can provide fundamental information for geophysics,geodynamics,seismology,and mineral exploration.Rectangular harmonic analysis(RHA)is proposed for regional gravity field modeling in this paper.By solving the Laplace’s equation of gravitational potential in local Cartesian coordinate system,the rectangular harmonic expansions of disturbing potential,gravity anomaly,gravity disturbance,geoid undulation and deflection of the vertical are derived,and so are the formula for signal degree variance and error degree variance of the rectangular harmonic coefficients(RHC).We also present the mathematical model and detailed algorithm for the solution of RHC using RHA from gravity observations.In order to reduce the edge effects caused by periodic continuation in RHA,we propose the strategy of extending the size of computation domain.The RHA-based modeling method is validated by conducting numerical experiments based on simulated ground and airborne gravity data that are generated from geopotential model EGM2008 and contaminated by Gauss white noise with standard deviation of 2 mGal.The accuracy of the 2.5′×2.5′geoid undulations computed from ground and airborne gravity data is 1 and 1.4cm,respectively.The standard error of the gravity disturbances that downward continued from the flight height of 4 km to the geoid is only 3.1 mGal.Numerical results confirm that RHA is able to provide a reliable and accurate regional gravity field model,which may be a new option for the representation of the fine structure of regional gravity field.     

等效偶层法位场曲面延拓的原理和计算方法

作者在本文中提出并证明Fredholm第二类积分方程幂级数解新的收敛条件。根据新的收敛条件,详细设计了等效偶层法位场曲面延拓的计算方法。在两种地形的试算中,证明了收敛条件的正确性和计算方法的有效性,对数值计算中的问题也进行了讨论。     

三维复杂山地多级次群推进迎风混合法多波型走时计算

三维复杂山地条件下的各种地震波型的走时计算技术,可以直接用于复杂山地区域地震波运动学特性的分析、地震数据采集观测系统的设计以及直接基于三维复杂地表的地震数据处理技术的研发.为了在三维复杂地表条件下准确、灵活且稳定地计算各种地震波型的走时,提出一种多级次群推进迎风混合法.该算法利用不等距迎风差分法简洁稳定地处理三维复杂地表及附近的局部走时计算问题,利用计算精度不错的迎风双线性插值法处理绝大部分均匀正方体网格中的局部走时计算问题,利用群推进法模拟三维复杂地表条件下地震波前的扩展问题,利用多级次算法处理各种类型的地震波的走时计算问题.算法分析和计算实例表明:新方法具有很好的计算精度与效率,且能灵活稳定地处理三维复杂地表复杂介质条件下的多波型走时计算问题.     

复杂地形条件下航空伽玛能谱地形改正方法探讨

应用矩形辐射体航空伽玛辐射场理论,研究应用于复杂地形条件下的航空伽玛能谱资料的地形改正方法,并针对该方法进行验证. 本文利用航空物探测量过程中获得的DTM数据（达到了地形改正所需的地形起伏数据精度）,应用矩形辐射体航空伽玛辐射场理论,根据地面辐射体与航空伽玛场分布之间的正演关系,对航空伽玛能谱解释方法的原理进行了正演分析和反演推导. 探索出一种按影响角进行地形改正的方法. 该方法特点是适合于任意飞行方式(缓地形和水平飞行均可)和任意地形条件航空伽玛能谱的逐点地形改正. 结果表明该地形改正方法能够基本消除航空伽玛能谱测量中的地形起伏产生的影响,经地形修正后的航空伽玛能谱异常能较正确地反映地面辐射体的真实情况.     

Three-dimensional topographic responses in MT using finite difference method

ＩｎｔｒｏｄｕｃｔｉｏｎＴｈｅｒｅａｒｅｍａｎｙｎｕｍｅｒｉｃａｌｍｏｄｅｌｉｎｇｓｔｕｄｉｅｓｏｆｔｏｐｏｇｒａｐｈｙｆｏｒＭＴｆｉｅｌｄｓａｐｐｅａｒｅｄｉｎｔｈｅｐａｓｔｌｉｔｅｒａｔｕｒｅ,ｂｕｔｍｏｓｔｏｆｔｈｅｗｏｒｋｓａｒｅ２－Ｄｒｅｓ...     

Effect of Near-fault Terrain upon Dislocation Modeling

—Coseismic surface deformation provides important information needed to determine source rupture geometry and slip distribution as well as to estimate seismic moment. In this study, numerical experiments were designed to analyze and classify how free-surface topography affects surface deformation. The investigation was performed by 3-D finite element modeling. Results of this study show that crustal deformation induced by near-fault terrain is significant and can be measured with present geodetic survey techniques. The characteristics of the terrain effects show that a hill structure produces more crustal deformation than a half-space model, and that the crustal deformation of a basin structure is less than that of the half-space model. The topographic correction is in the order of five percent of the fault dislocation. On the basis of the relationship between fault offset and earthquake magnitude, it is suggested that the terrain effects on the coseismic crustal deformation of shallow earthquakes with a magnitude greater than 5.6 should be considered as one of the major errors in coseismic deformation modeling which ignored the surface topography on the order of 300 meters.     

A robust channel network extraction method combining discrete curve evolution and the skeleton construction technique

The automatic mapping of drainage networks from terrain representation has been an interesting topic in hydrological and geomorphological modeling. However, the existing methods often suffer from high sensitivity to terrain noise or lose significant stream branches and accurate channel paths. In this paper, we propose a contour-based framework in drainage network extraction. The proposed framework incorporates discrete curve evolution (DCE) to eliminate the noise influence by dynamically segmenting the contour lines (CLs) into valley bends, and to detect the valley feature points. The skeleton construction technique is then applied to distill more accurate channel paths in complex terrain. Finally, a linking step is undertaken to generate the channel network. The proposed method was tested on a series of elevation datasets, with varied resolution, region size, and local relief. The experiments verified that the proposed method can achieve highly accurate channel networks and is robust, even in regions with high-contrast relief, and/or in cases with significant terrain noise and irregularities.