Numerical computation of spherical harmonics of arbitrary degree and order by extending exponent of floating point numbers

By extending the exponent of floating point numbers with an additional integer as the power index of a large radix, we compute fully normalized associated Legendre functions (ALF) by recursion without underflow problem. The new method enables us to evaluate ALFs of extremely high degree as 2³² =  4,294,967,296, which corresponds to around 1 cm resolution on the Earth’s surface. By limiting the application of exponent extension to a few working variables in the recursion, choosing a suitable large power of 2 as the radix, and embedding the contents of the basic arithmetic procedure of floating point numbers with the exponent extension directly in the program computing the recurrence formulas, we achieve the evaluation of ALFs in the double-precision environment at the cost of around 10% increase in computational time per single ALF. This formulation realizes meaningful execution of the spherical harmonic synthesis and/or analysis of arbitrary degree and order.     

A unified approach to the Clenshaw summation and the recursive computation of very high degree and order normalised associated Legendre functions

 Spherical harmonic expansions form partial sums of fully normalised associated Legendre functions (ALFs). However, when evaluated increasingly close to the poles, the ultra-high degree and order (e.g. 2700) ALFs range over thousands of orders of magnitude. This causes existing recursion techniques for computing values of individual ALFs and their derivatives to fail. A common solution in geodesy is to evaluate these expansions using Clenshaw's method, which does not compute individual ALFs or their derivatives. Straightforward numerical principles govern the stability of this technique. Elementary algebra is employed to illustrate how these principles are implemented in Clenshaw's method. It is also demonstrated how existing recursion algorithms for computing ALFs and their first derivatives are easily modified to incorporate these same numerical principles. These modified recursions yield scaled ALFs and first derivatives, which can then be combined using Horner's scheme to compute partial sums, complete to degree and order 2700, for all latitudes (except at the poles for first derivatives). This exceeds any previously published result. Numerical tests suggest that this new approach is at least as precise and efficient as Clenshaw's method. However, the principal strength of the new techniques lies in their simplicity of formulation and implementation, since this quality should simplify the task of extending the approach to other uses, such as spherical harmonic analysis. Received: 30 June 2000 / Accepted: 12 June 2001     

Geometric Algebra Model for Geometry‐oriented Topological Relation Computation

Classical topological relation expressions and computations are primarily based on abstract algebra. In this article, the representation and computation of geometry‐oriented topological relations (GOTR) are developed. GOTR is the integration of geometry and topology. The geometries are represented by blades, which contain both algebraic expressions and construction structures of the geometries in the conformal geometric algebra space. With the meet, inner, and outer products, two topology operators, the MeetOp and BoundOp operators, are developed to reveal the disjoint/intersection and inside/on‐surface/outside relations, respectively. A theoretical framework is then formulated to compute the topological relations between any pair of elementary geometries using the two operators. A multidimensional, unified and geometry‐oriented algorithm is developed to compute topological relations between geometries. With this framework, the internal results of the topological relations computation are geometries. The topological relations can be illustrated with clear geometric meanings; at the same time, it can also be modified and updated parametrically. Case studies evaluating the topological relations between 3D objects are performed. The result suggests that our model can express and compute the topological relations between objects in a symbolic and geometry‐oriented way. The method can also support topological relation series computation between objects with location or shape changes.     

Global spherical harmonic computation by two-dimensional Fourier methods

A method is presented for performing global spherical harmonic computation by two-dimensional Fourier transformations. The method goes back to old literature (Schuster 1902) and tackles the problem of non-orthogonality of Legendre-functions, when discretized on an equi-angular grid. Both analysis and synthesis relations are presented, which link the spherical harmonic spectrum to a two-dimensional Fourier spectrum. As an alternative, certain functions of co-latitude are introduced, which are orthogonal to discretized Legendre functions. Several independent Fourier approaches for spherical harmonic computation fit into our general scheme.     

多种超高阶次缔合勒让德函数计算方法的比较

扩展高阶和超高阶重力场模型的构制与应用的数值稳定性取决于超高阶次缔合勒让德函数的计算方法。文中详细介绍了现有的多种缔合勒让德函数的递推计算方法:标准前向列推法、标准前向行推法、跨阶次递推法和Belikov列推法。从计算速度、计算精度和计算溢出问题3个角度分析比较了阶次高至2 160阶的各种方法的优劣。通过数值试验证明,Belikov列推法和跨阶次递推法是计算超高阶次缔合勒让德函数较优的方法,而其他几种方法不能用于超高阶次缔合勒让德函数的计算。文中结论为超高阶次球谐综合与球谐分析的数值计算提供了可靠的依据。     

Truncation error formulae for the disturbing gravity vector

Recurrence relations have been derived for truncation error coefficients of the extended Stokes' function and its partial derivatives required in the computation of the disturbing gravity vector at any elevation above the earth's surface. The corresponding formulae, the example of values of the truncation error coefficients for H=30.1 km and ψ₀=30^∘ and the estimations of truncation error are given in this article. Received: 26 January 1996 / Accepted: 11 June 1997     

Recurrence relations for the truncation error coefficients for the extended stokes function

Recurrence relations for the truncation error coefficients of the extended Stokes function required in the computation of gravimetric geoidal heights at any elevation above the earth's surface have been derived. The computation of these coefficients generally involves a small fixed number of terms except at altitudes of2700 km or more when one of the terms involved has to be computed from an infinite series. To confirm the accuracy of the coefficients a verification formula has been devised which uses a series expansion of a piece-wise continuous function such that it is equal to the extended Stokes function over a given range but vanishes elsewhere. Contribution of the Earth Physics Branch #1053.     

A new formula for evaluating the truncation error coefficient

In this paper, a new formula for evaluating the truncation coefficientQ _n is derived from recurrence relations of Legendre polynomials. The present formula has been conveniently processed by an electronic computer, providing the value ofQ _n up to a degreen=49 which are exactly equal to those of Paul (1973).     

The Generalized Conic Projection

The Equidistant Cylindric and the Equidistant Conic projections are reviewed as foundations for a series of new projections. The new projections are introduced as the Patchwork Conic, an amalgamation of Equidistant Conic projections; the Generalized Equidistant Cylindric, with a standard path, or conformal path of correct scale, that is definable by the projection designer; and the Generalized Equidistant Conic, having two definable standard paths. Several potential uses are mentioned. Maps, formulae, and computation examples are provided. It is concluded that the projection is unique in allowing tailoring without compromising desirable mathematical properties.     

Formation of surface spherical harmonic normal matrices and application to high-degree geopotential modeling

 The structure of normal matrices occurring in the problem of weighted least-squares spherical harmonic analysis of measurements scattered on a sphere with random noises is investigated. Efficient algorithms for the formation of the normal matrices are derived using fundamental relations inherent to the products of two surface spherical harmonic functions. The whole elements of a normal matrix complete to spherical harmonic degree L are recursively obtained from its first row or first column extended to degree 2L with only O(L ⁴) computational operations. Applications of the algorithms to the formation of surface normal matrices from geoid undulations and surface gravity anomalies are discussed in connection with the high-degree geopotential modeling. Received: 22 March 1999 / Accepted: 23 December 1999     

Sign-constrained robust least squares, subjective breakdown point and the effect of weights of observations on robustness

The findings of this paper are summarized as follows: (1) We propose a sign-constrained robust estimation method, which can tolerate 50% of data contamination and meanwhile achieve high, least-squares-comparable efficiency. Since the objective function is identical with least squares, the method may also be called sign-constrained robust least squares. An iterative version of the method has been implemented and shown to be capable of resisting against more than 50% of contamination. As a by-product, a robust estimate of scale parameter can also be obtained. Unlike the least median of squares method and repeated medians, which use a least possible number of data to derive the solution, the sign-constrained robust least squares method attempts to employ a maximum possible number of good data to derive the robust solution, and thus will not be affected by partial near multi-collinearity among part of the data or if some of the data are clustered together; (2) although M-estimates have been reported to have a breakdown point of 1/(t+1), we have shown that the weights of observations can readily deteriorate such results and bring the breakdown point of M-estimates of Huber’s type to zero. The same zero breakdown point of the L ₁-norm method is also derived, again due to the weights of observations; (3) by assuming a prior distribution for the signs of outliers, we have developed the concept of subjective breakdown point, which may be thought of as an extension of stochastic breakdown by Donoho and Huber but can be important in explaining real-life problems in Earth Sciences and image reconstruction; and finally, (4) We have shown that the least median of squares method can still break down with a single outlier, even if no highly concentrated good data nor highly concentrated outliers exist. An erratum to this article is available at .     

不同扰动位泛函间的积分变换广义核函数

基于物理大地测量边值问题的解,利用一阶边界算子定义,推导重力异常Δg、单层密度μ、大地水准面高N,垂线偏差ε、扰动重力δg等扰动场元的解。利用球谐函数的正交特性,通过对核函数的算子运算,可以得到上述扰动场元的有关逆变换公式。相对经典物理大地测量公式应用的边界面条件,笔者将含有因子r的对应扰动场元反演关系的公式称为广义积分公式。针对常用的重力异常Δg、大地水准面高N,垂线偏差ε、扰动重力δg计算,重点分析它们之间的变换关系,给出利用某个选定扰动场元计算其他扰动场元的广义积分公式。同时,通过对积分边界面的讨论,分析经典公式与广义积分公式的差异和联系。最后,给出所有外部扰动场元与核函数映射的关系表。     

基于V9I的空间关系映射与操作

利用V9I模型中目标的边界、内部和Voronoi区域均可量测及易于操作的特点 ,研究建立底层数据结构与空间关系语义层之间的V9I映射机制 ,在Voronoi的动态栅格生成算法的基础上 ,构建空间关系的基本操作。最后 ,以VC 为开发工具并采用面向对象的技术 ,设计了基于V9I的空间关系操作工具原型VTKit ,并给出部分操作实例     

Topological relations between spherical spatial regions with holes

ABSTRACT

There is growing interest in globally modelling the entire planet. Although topological relations between spherical simple regions and topological relations between regions with holes in the plane have been investigated, few studies have focused on the topological relations between spherical spatial regions with holes. The 16-intersection model (16IM) is proposed to describe the topological relations between spatial regions with holes. A total of 25 negative conditions are proposed to eliminate the impossible topological relations between spherical spatial regions with holes. The results show that (1) 3 disjoint relations, 3 meet relations, 66 overlap relations, 7 cover relations, 3 contain relations, 1 equal relation, 7 coveredBy relations, 3 inside relations, 1 attach relation, 52 entwined relations, and 28 embrace relations can be distinguished by the 16IM and that (2) the formalisms of attach, entwined, and embrace relations between the spherical spatial regions without holes based on the 9IM and that between the spherical spatial regions with holes based on the simplified 16IM are different, whereas the formalisms of other types of relations between spherical spatial regions without holes based on the 9IM and that between the spherical spatial regions with holes based on a simplified 16IM are the same.     

LEAST-SQUARE SOLUTIONS WITH WEIGHTS

Abstract

The modern abridged method of solution as applied to observation equations was given in an earlier number of this Review. The present article, applying the same abridged method to the solution of conditioned equations, shows how the weights of any functions of the adjusted values can be obtained as well as the corrections themselves. The case where the weights of adjusted functions are most frequently required in practice is that of a triangulation base extension figure. It is well known that the error generated in reaching the first side of the main triangulation proper from the measured base may be considerably greater than thereafter when sides are long and grazing shots rare. For this reason it is good practice to compute rigorously the probable error of the first main side from the base. It is usually found that the p.e. of the base itself is only a small fraction of the total. If the total p.e. of this first main side is too large, it may then be considered whether the extension figure should not be either completely reobserved or even redesigned to give a better result. Even if the observations themselves are all considered to be of equal weight (as is common practice nowadays with good instruments and good methods of observing), the weights of the adjusted functions will still differ from the mean, depending on the condition equations set up by the extension figure adopted.     

On computing ellipsoidal harmonics using Jekeli’s renormalization

Gravity data observed on or reduced to the ellipsoid are preferably represented using ellipsoidal harmonics instead of spherical harmonics. Ellipsoidal harmonics, however, are difficult to use in practice because the computation of the associated Legendre functions of the second kind that occur in the ellipsoidal harmonic expansions is not straightforward. Jekeli’s renormalization simplifies the computation of the associated Legendre functions. We extended the direct computation of these functions—as well as that of their ratio—up to the second derivatives and minimized the number of required recurrences by a suitable hypergeometric transformation. Compared with the original Jekeli’s renormalization the associated Legendre differential equation is fulfilled up to much higher degrees and orders for our optimized recurrences. The derived functions were tested by comparing functionals of the gravitational potential computed with both ellipsoidal and spherical harmonic syntheses. As an input, the high resolution global gravity field model EGM2008 was used. The relative agreement we found between the results of ellipsoidal and spherical syntheses is 10^?14, 10^?12 and 10^?8 for the potential and its first and second derivatives, respectively. Using the original renormalization, this agreement is 10^?12, 10^?8 and 10^?5, respectively. In addition, our optimized recurrences require less computation time as the number of required terms for the hypergeometric functions is less.     

New characteristics of weighted GDOP in multi-GNSS positioning

In positioning, navigation and timing applications of multi-GNSS (global navigation satellite system) constellations, the geometric dilution of precision (GDOP) offers an important index for selecting satellites and evaluating positioning accuracy. However, GDOP assumes that the measurement errors of all the tracked satellites are independent and have the same accuracy level, which is impossible in practice, especially when the tracked satellites are from various constellations. Through introducing a weighted matrix describing the measurement errors of different satellites into a common GDOP, we focus on new characteristics of weighted GDOP (WGDOP) in two aspects. First, we compare the sizes of WGDOP and the common GDOP based on the range of the weights of different satellites, i.e., the diagonal elements of the weighted matrix. In addition, when the weights of different satellites increase, the change of WGDOP with the weights is also derived. Moreover, a closed-form formula for calculating WGDOP is also presented. The theoretical derivations demonstrate that the closed-form can reduce the computation burden effectively. Furthermore, numerical tests verify these analyses.     

DIRECTIONS AND ANGLES IN A SIMPLE CHAIN

Abstract

Mr. Rainsford's article on “Least Square Adjustments of Triangulation: Directions versus Angles” in the Empire Survey Review No. 78, Vol. x, October 1950, leads to many speculations and interesting results. I try to show here, how, by assuming artifices to simplify the results, weights may be assigned to angles derived from directions so that the results of adjustment by angles, with these weights, will be the same as the adjustment by directions, all of equal weight.     

Recurrence relations for integrals of Associated Legendre functions

Recurrence relations for the evaluation of the integrals of associated Legendre functions over an arbitrary interval within (0°, 90°) have been derived which yield sufficiently accurate results throughout the entire range of their possible applications. These recurrence relations have been used to compute integrals up to degree 100 and similar computations can be carried out without any difficulty up to a degree as high as the memory in a computer permits. The computed values have been tested with independent check formulae, also derived in this work; the corresponding relative errors never exceed 10⁻²³ in magnitude. Contribution from the Earth Physics Branch No. 719     

Strength analysis of geodetic control networks

In strength analysis of horizontal geodetic networks it is appropriate to use pairs of functions which involve the relative position of two points and relative position of three points. Using properly chosen pairs of functions, formulae are given which allow the computation of precision criteria for the orientation and scale of the network as well as its shape. To illustrate the presentation of results, new types of errors ellipses are introduced. Analogies and correlations existing among the adopted functions are introduced by using the concept of orthogonal networks which are defined in the paper.