Full text
Journal of Geodesy (2024) 98:47 https://doi.org/10.1007/s00190-024-01854-1 ORIGINAL ARTICLE A novel geometric method based on conformal geometric algebra applied to the resection problem in two and three dimensions Jorge Ventura1·Fernando Martinez1·Francisco Manzano-Agugliaro1·Aleš Návrat2·Jaroslav Hrdina2· Ahmad H. Eid3·Francisco G. Montoya1 Received: 5 October 2023 / Accepted: 23 April 2024 / Published online: 27 May 2024 © The Author(s) 2024 Abstract This paper introduces a novel method for solving the resection problem in two and three dimensions based on conformal geometric algebra (CGA). Advantage is taken because of the characteristics of CGA, which enables the representation of points, lines, planes, and volumes in a unified mathematical framework and offers a more intuitive and geometric understanding of the problem, in contrast to existing purely algebraic methods. Several numerical examples are presented to demonstrate the efficacy of the proposed method and to compare its validity with established techniques in the field. Numerical simulations indicate that our vector geometric algebra implementation is faster than the best-known algorithms to date, suggesting that the proposed GA-based methods can provide a more efficient and comprehensible solution to the twoand three-dimensional resection problem, paving the way for further applications and advances in geodesy research. Furthermore, the method’s emphasis on graphical and geometric representation makes it particularly suitable for educational purposes, allowing the reader to grasp the concepts and principles of resection more effectively. The proposed method has potential applications in a wide range of other fields, including surveying, robotics, computer vision, or navigation. Keywords Resection problem ·Triangulation ·Snellius–Pothenot ·Conformal geometric algebra BFrancisco Manzano-Agugliaro [email protected] Jorge Ventura [email protected] Fernando Martinez [email protected] Aleš Návrat na[email protected].cz Jaroslav Hrdina [email protected].cz Ahmad H. Eid [email protected] Francisco G. Montoya [email protected] 1Department of Engineering, University of Almeria, Crtra. Sacramento s/n, 04120 Almería, Almería, Spain 2Institute of Mathematics, Faculty of Mechanical Engineering, Brno University of Technology, Technická 2896/2, Královo Pole, 61669 Brno, Czech Republic 3Electrical Engineering Department, Port Said University, Port-Said, Egypt 1 Introduction The resection problem, also known in surveying as the Snellius–Pothenot (SP) or the inverse intersection problem, involves calculating the position of an unknown point P(also called a station) using the positions of three known points A, B, and C, and relative angular measurements from P.Itisa relevant problem not only in geodesy and surveying, but also in other disciplines such as robot path planning (Masselli and Zell 2014), positioning (Pierlot and Van Droogenbroeck 2014), navigation (Pereira et al. 2018) or computer graphics (Mazaheri and Habib 2015), and can be solved both geometrically and algebraically. The solution of the two-dimensional problem has been known for several centuries and has numerous variants [almost 500! according to Bock (1956)]. The three-dimensional problem is much more intricate and only complex and sophisticated algebraic solutions are known. The twoand three-dimensional configurations are illustrated in Fig. 1. 123
47 Page 2 of 21 J. Ventura et al. Fig. 1 Representation of the resection problem in 2D (left) and 3D (right) 1.1 Motivation The determination of the position of an observer based on angular measurements from known points is of interest in several disciplines, such as surveying, computer graphics, optics, and robotics. Traditionally, solutions to this problem have relied on heavily algebraically loaded methods, which can be complex and challenging to comprehend. Furthermore, these methods do not always provide an intuitive understanding of the geometric relationships involved in the problem. In both two-dimensional (2D) and three-dimensional (3D) problems, there is a need for a more geometric or graphical approach that simplifies the study of the resection problem. Although there are existing graphical methods for solving the problem in 2D that have been known for some time, their algebraic implementation can be quite cumbersome, hindering their widespread adoption. When the problem is approached from a geometric perspective, a better understanding of the underlying structures and relationships can be achieved, making the problem more accessible to a wider range of researchers and practitioners. This could potentially lead to novel applications and advancements in the related areas mentioned above. In light of the exposed ideas, the main motivation behind this paper is to develop a novel geometric method based on conformal geometric algebra (CGA) to address the resection problem in twoand three-dimensional ways. By leveraging the unifying properties of CGA to represent geometric primitives within a single mathematical framework, we aim to provide a more intuitive and geometric understanding of the problem, simplifying its study, and paving the way for further applications and advancements in the field. 1.2 Background and literature overview The resection problem has been extensively studied in the literature, with various algebraic and geometric or graphical methods proposed for its solution. It appears that ancient Greeks, such as Hipparchus or Ptolemy, already studied this problem in the context of astronomy, although the first person to solve the problem in the context of surveying was the Dutch mathematician Willebrord Snel van Royen (known Table 1 Classification of main methods to solve the 2D resection problem (see Pierlot and Van Droogenbroeck 2014) Method Group Year Snellius Trigonometric 1617 Collins Geometric/Graphical 1671 Kaestner-Burkhardt Trigonometric 1801 McGillem Geometric 1988 Cohen and Koss Geometric 1992 Madsen and Andersen Trigonometric 1998 Sanchiz et al Iterative 2004 Esteves et al Geometric 2006 Easton and Cameron Trigonometric 2006 Tsukiyama Geometric 2009 Font-Llagunes Geometric 2009 Dekov Iterative 2012 Ligas Geometric 2013 Pierlot et al Geometric 2014 Willerding Other 2020 Tienstra Other Unknown Cassini Geometric/Trigonometric Unknown as Snellius) in 1617 (Wreede 2007). He achieved the goal by using geometric and algebraic methods mainly based on trigonometry. This same problem was addressed in 1671 by John Collins in his work Philosophical Transactions. Collins contributed significantly to the topic by presenting a new and elegant geometrical solution, which involves the use of just one circle and an auxiliary point. In 1692, Laurent Pothenot, who was working on the definition of the meridian north of Paris, presented a paper on the subject. However, according to McCaw (1918) and others, Pothenot did not contribute anything new to the solution of the problem and all he did was publish the works of Snellius and Collins under his own name. Other authors who have studied this issue are shown in Table 1. In the context of geodesy and surveying, the resection problem has been a fundamental problem for centuries. The problem has been addressed from various perspectives, including algebraic (Awange et al. 2010), geometric/graphical (Masselli and Zell 2014), and numerical 123
A novel geometric method based on conformal geometric algebra applied to the resection problem… Page 3 of 21 47 (Dekov 2012) methods. Algebraic solutions to the problem are well established and have been extensively studied. However, these methods often involve complex algebraic manipulations and do not provide any intuitive geometric understanding of the problem. The two-dimensional problem has been solved using approaches such as graphical methods (Bock 1956), analytical geometry (Bil 1992), matrix methods (Bähr 1991), and algebraic methods based on Sylvester resultants (Awange and Grafarend 2002), Groebner bases (Awange 2002), and other elimination techniques (Sturmfels 2002). Awange and Grafarend (2005) provided a detailed overview of several algebraic techniques for solving 2D and 3D resection. The three-dimensional version is more complex, with most solutions relying on intricate algebraic methods. Grunert originally formulated the problem in 1841 and derived quartic equations to determine the unknown distances (Grunert 1841). Since then, numerous procedures have been developed to optimise Grunert’s formulation and reduce the computational steps involved (Awange and Grafarend 2005; Fischler and Bolles 1981). Algebraic techniques such as Groebner bases (Awange and Grafarend 2003a), polynomial resultants (Awange and Grafarend 2003b), and linear homotopy (Paláncz et al. 2010) have been effectively applied to solve the 3D problem. However, geometric solutions to the resection problem have been less explored. The geometric or graphical approach provides a more intuitive understanding of the problem and can be more easily visualised. In recent years, there has been a growing interest in the application of geometric algebra (GA) to solve several geometric problems. GA provides a unified mathematical framework for representing and manipulating geometric objects, making it a powerful tool for geometric computations (Dorst et al. 2010). The application of GA to the resection problem was first proposed by Smith (2023b,a). However, the application of GA to the resection problem is still in its early stages, and there is still a long way to go. In particular, the application of conformal geometric algebra (CGA) to the resection problem has not been fully explored. CGA extends GA by incorporating the concept of conformal transformations, which can provide a more powerful and flexible framework for geometric computations (Hestenes and Sobczyk 2012; Dorst et al. 2010). CGA provides a unified framework for handling points, lines, planes, circles, spheres, and other geometric entities. Each entity has a unique representation, and the relationships between entities correspond to algebraic relationships between their CGA representations. This eliminates the need for coordinate transformations when transitioning between reference frames. Cameron and Lasenby (2008) showed that CGA subsumes projective geometry and is more computationally efficient than matrix methods. CGA has found applications in computer vision (Wareham et al. 2004), robotics (Zaplana et al. 2022), and other geometric computing problems (Hitzer et al. 2022). However, its potential for solving resection-type problems in geodesy and surveying remains relatively unexplored. This paper seeks to fill the void by introducing a new geometric solution based on CGA for the 2D and 3D resection problem. 1.3 Contributions In this paper, a new method for solving the resection problem using the mathematical framework geometric algebra in 2D and 3D dimensions is proposed. The method provides a simple solution based on purely geometric and graphical principles. Specifically, the 3D version of the problem is thoroughly analysed in detail. Specific novel contributions are the following: •A novel conformal geometric algebra (CGA)-based method is presented to solve the twoand threedimensional resection problem, providing a more intuitive and unified geometric approach. •The proposed method is compared with established techniques, and its advantages and efficacy are demonstrated through numerical examples. •The potential applications of the method in several fields such as computer graphics, optics, and robotics are highlighted, emphasising its versatility, and paving the way for future advancements in geometric research. •Several algorithms have been developed that provide better results than the best-known algorithms to date from a computational perspective. 1.4 Outline The remainder of the article is organised as follows. Section 2 provides an introduction to geometric algebra and conformal geometric algebra. Section3revisits the resection problem and reviews traditional methods. Section 4presents existing and new GA-based methods. Section 4presents the proposed method based on CGA to solve the resection problem. Section5provides several application examples to demonstrate the effectiveness of the proposed method, along with benchmarks for computational efficiency. Section 6provides an error and uncertainty analysis of the proposed methods. Finally, Sect. 7concludes the article with a summary and some suggestions for future work. 2 Basics concepts in geometric algebra Geometric algebra is a mathematical framework for representing geometric objects and transformations in a unified 123
47 Page 4 of 21 J. Ventura et al. way (Hestenes and Sobczyk 2012). GA extends the algebra of vectors to include other geometric objects such as points, lines, planes, and volumes. GA provides a powerful tool to solve geometric problems and has applications in a wide range of fields, including electrical engineering (Montoya et al. 2019,2021), computer vision (Hrdina and Návrat 2017), robotics (Hrdina et al. 2017), and other engineering fields [see Hitzer et al. (2022) and references there in]. GAs are used mainly in situations where Euclidean transformations play a significant role. A simple way to introduce GAs is to understand them through more familiar tools such as complex numbers or quaternions (indeed, they are subalgebras of GA). For example, implementing Euclidean rotations using quaternions is essential to improve computational capabilities and adopt an object-oriented approach. We can see quaternions Has a natural extension of the complex numbers Cin the form z=a+bi+cj+dk(1) where a,b,c,d∈Rand the products of the basis elements i,j,ksatisfy the multiplications rules i2=j2=k2=−1 and ij =−ji =k. The quaternion Im(z)=bi+cj+dk is called the imaginary part of z, and Re(z)=ais called the real part of z. Note the use of bold letters for quaternions and regular font for real numbers. The operation of rotating an object by an angle θaround the axis n=nxi+nyj+nzkcan be represented by the quaternion R=eθ 2(nxi+nyj+nzk) =1 2cos θ+nxi+nyj+nzksin θ(2) acting on the vector q=xi+yj+zk. Note that ||n|| = √n¯ n=1, where the bar decoration stands for quaternionic conjugation. The representation is done in a similar way as in the case of complex numbers. The difference is that quaternions are not commutative and thus they act via the so-called sandwich product RqR−1=1 ||R||Rq ¯ R(3) So, we are working in a four-dimensional linear space. If we want to realise also the translations, we have to extend the algebra by another dimension using an element such as 2=0 and i=j=k=0. This algebra is usually called dual quaternions. For more on the use of quaternions in general engineering topics, see Selig (2005). On the other hand, the wedge operation on a vector space allows us to work with linear subspaces. A line can be characterised by a vector ,soxbelongs to a line if and only if x∧=0. The wedge of two vectors then characterises the plane in the same way. A vector space closed in the wedge operation is called a Grassmannian algebra. The combination of these two concepts leads to the notion of GA. 2.1 Euclidean vector GA Vector geometric algebra (VGA) is perhaps one of the simplest GAs. For the two-dimensional case, the VGA (G2)isa Grassmannian algebra based on two orthonormal basis vectors (σ1,σ2) together with a bilinear operation known as geometric product satisfying the following identities: σ1σ2=σ1·σ2+σ1∧σ2=σ1∧σ2 σ1·σ2=0 σ2 1=σ2 2=1(4) The term σ1∧σ2is known as bivector (σ12 for short). Geometrically, a bivector represents an oriented plane segment spanned by the two vectors. We can use bivectors to represent rotations in the following way. The collection of bivectors σ1∧σ2forms a one-dimensional vector space that is closed under multiplication. We can then generate rotations by applying the exponential map to bivectors. For example, the exponential of the bivector e1 2θσ12 =1 2(cos θ+sin θσ12)(5) is used as in (3) to perform a rotation by an angle θin the plane spanned by σ1and σ2. In this way, the group of rotations generated by the bivectors is mathematically equivalent to the group of unitary complex numbers eiθ, which also represent rotations in the complex plane (note that (σ1∧σ2)2=−1). However, bivectors provide a more intuitive geometric representation of rotations directly in vector space. Similarly, we can introduce the G3algebra for the purposes of reasoning in the 3D space. The algebra G3is based on three generators σ1,σ2and σ3together with a geometric product defined by the following identities: σiσj=σi∧σjwhere i= j σ2 i=1 where i=1,2,3(6) In Sect.4.1, it will be shown how to use the VGA-based method to solve a 2D version of the resection problem. This procedure is, in fact, a use of GA G2. 2.2 Conformal geometric algebra The goal here is to create a model of Euclidean geometry. Specifically, geometry whose symmetry group contains the 123
A novel geometric method based on conformal geometric algebra applied to the resection problem… Page 5 of 21 47 Euclidean symmetries (rotations, translations, etc.). For that purpose, a nondegenerate quadratic form will be chosen, thus obtaining the conformal geometric algebra (CGA). In the case of CGA for a two-dimensional Euclidean space, a GA of signature (3,1)is obtained [also known as Compass Rule Algebra CRA, see Hildenbrand (2018)], with basis vectors σ1,σ2,σ+and σ−such that σ2 i=1,i∈{1,2,+},σ2 −=−1(7) σiσj=−σjσi,i= j,i,j∈{1,2,+,−} (8) For mathematical convenience, it is advisable to define two new basis vectors σ0=1 2(σ−−σ+)(9) σ∞=σ−+σ+,(10) with properties σ2 0=σ2 ∞=0 and σ0·σ∞=−1 (11) Vector σ0is known to represent the Euclidean point at the origin of the coordinate system, and vector σ∞represents the point at infinity. Using the rules above, it can be proved that any Euclidean point Xcan be mapped to the CGA vector space as X−→ x=σ0+(x1σ1+x2σ2)+1 2(x2 1+x2 2)σ∞ (12) In “Appendix A.1”, explicit calculations that justify our choices can be found. As in Sect.2.1, the extension to higher dimensions is straightforward by adding the element σ3.This will make G3appear instead of G2, and quaternions will appear in the bivectors. In the notation used throughout this paper to ensure clarity, regular font indicates real numbers, boldface denotes vectors in the CGA space, boldface with an overhead arrow (e.g. x) represents vectors in the Euclidean space, and uppercase letters in boldface indicate multivectors in CGA. The usefulness of CGA and CRA is based on the fact that the distance between Euclidean points is encoded in the scalar product. The points Xand Yare represented by the Euclidean vectors xand y, respectively. The mapping defined in (12) maps these Euclidean points (vectors) to vector elements x and yin the CGA space. As shown in “Appendix A.2”, the scalar product between two CGA points codes the distance: x·y=−1 2 x− y2(13) Thus, the Euclidean point X(represented by the Euclidean vector x) lies on the sphere (circle) Swith centre in point C and radius rif and only if it satisfies the identity x·c=−1 2r2 in the conformal space. Since x·σ∞=−1 for each point, this identity may be written as x·c+1 2r2=x·c−1 2r2(x·σ∞) =x·c−1 2r2σ∞=0(14) Thus, the vector element c−1 2r2σ∞represents the sphere Sin the CGA space. It should be emphasised that within CGA, there exist two distinctive methods to reference the identical geometric entity: IPNS and OPNS (see “Appendix A.3”). IPNS is superior for transformations and intersections, whereas OPNS is advantageous for blending and morphing tasks. It is crucial to be aware of the representation we are operating in; however, transitioning between representations can be achieved seamlessly using the dual operator (), as described in “Appendix A.3”. 3 Revisiting the resection problem In the context of geodesy and surveying, the resection problem plays an important role in determining the position of an observer based on angular measurements from three known reference points. Over the years, various approaches have been developed to address this problem, ranging from graphical or geometrical methods to algebraic ones. However, many of these methods can be complex or tedious, particularly when addressing the problem in three dimensions. With the growing importance of accurate positioning in modern applications, it is essential to revisit the resection problem and explore innovative approaches that offer more intuitive solutions. 3.1 Traditional methods In this section, traditional methods for solving the resection problem are reviewed. They can be classified into four basic groups: trigonometric, geometric/graphical, iterative (numerical) and others (see Table 2). A comprehensive list of methods for 2D resection problems is already presented in Table 1. Note that sometimes the frontier between trigonometric and geometric methods is not as clear because both solutions exist at the same time (e.g. the Cassini method). Trigonometric methods are some of the oldest and most famous procedures for solving the 2D resection problem. They are based on the use of trigonometric functions to compute the position of the observer using the angles between the known points and the observer. 123
47 Page 6 of 21 J. Ventura et al. Table 2 Taxonomy of Resection Methods Group Methods Trigonometric Use trigonometric functions and equations to calculate station/robot position Geometric Compute intersection of circles or lines passing through known points and station/robot position Iterative Start with estimate of robot position, iteratively refine using error minimisation Others Use of Barycentric Coordinates or Complex Numbers The Pothenot–Snellius method, also known as the Kästner– Burkhardt method, is one of the most known and oldest procedures in this group. The Cassini method is another example of a trigonometric solution (with a graphical solution, too), which is similar to the one described by Esteves et al. and Cohen and Koss. The Easton and Cameron method is also a trigonometric approach, similar to the Cassini method. Geometric or graphical methods, on the other hand, use geometrical constructions and properties of the known points to calculate the position of the observer. These methods are based on the use of geometric and graphic principles to solve the problem. Esteves et al. proposed a geometric method that uses the intersection of circles to calculate the position of the observer. Cohen and Koss also proposed a geometric method that uses the intersection of circles, but with a different approach from Esteves et al. Font-Llagunes proposed a method that uses the intersection of lines and similarity between triangles to calculate the position of the observer. Pierlot et al. and Ligas proposed a method that uses the intersection of power lines to calculate the position of the observer, although the solution is given in algebraic form. Iterative methods use iterative algorithms to converge to the observer position. These methods are based on the use of an initial estimate of the observer position, which is refined iteratively until convergence is achieved. Iterative search is an example of an iterative method that uses a search algorithm to find the robot position. Sanchiz et al. proposed a method that uses an iterative search algorithm to find the observer position. Finally, there are other methods that do not fit the previous groups. The Tienstra method is one such example, which is a completely different approach based on barycentric coordinates. Another method like Willerding is based on the use of complex numbers to compute rotations in the Argand plane to find the observer position. 4 Resection using geometric algebra Nowadays, methods based on geometric algebra (GA) have been developed, providing a new solution to the resection problem while maintaining the focus on its geometrical roots. GA offers a versatile framework that can be adapted to different contexts based on the selection of specific metrics and the number of dimensions. For example, when all elements of the base σisquare to +1, the resulting algebra is known as vector geometric algebra (VGA) (see Sect.2.1). By extending the ability of the basis elements to square to −1 or 0 or by accommodating more dimensions, it becomes possible to explore alternative forms of GA. One such example is CGA, which incorporates two additional dimensions: one dimension squaring to +1 and another squaring to −1. This flexibility enables GA to solve a wide array of applications and problem domains, seamlessly scaling the number of dimensions in a straightforward way. 4.1 Vector GA method The 2D VGA-based method has been recently proposed by Smith (2023a) and published as disseminative material. The process is mainly geometric and results in obtaining a vector pthat describes the position of the point Pwhen choosing the middle point (B) as the origin. In this case, we start with a vector basis consisting of two elements σ={σ1,σ2}. Figure 2illustrates a representation of the problem, as well as a detailed sequence of steps carried out. First, using the known data (A,B,C,α, and β), circles c1and c2are drawn using points A,B,Pand B,C,P, respectively (see Fig. 2a). These circles serve as an auxiliary element for better understanding the solution, but are not required as such. With the help of the central angle theorem, the vectors d1and d2are obtained d1=v1+v1 tan α σ12 =v1 sin αe(90−α)σ12 d2=v2−v2 tan β σ12 =v2 sin βe(β−90)σ12 (15) where v1=A−Band v2=C−B. Note that Bwas chosen as the origin, but any other point can also be selected under the condition that αor βis not null. This situation occurs when Pis collinear with two of the three known points. In such a case, other points can be selected as the origin. Equation (15) indicates that vectors diare the result of rotating and scaling vectors vi, as shown in Fig. 2b. Note that the rotation angle is given by (90 −α) and (90 −β), respectively. The next step involves determining the vector das d2−d1. Finally, 123
A novel geometric method based on conformal geometric algebra applied to the resection problem… Page 7 of 21 47 Fig. 2 Vector GA-based method steps to solve the 2D resection problem (a) Initial setup, with (unknown) circles c1(A, Band P)andc2(B,Cand P). (b) Using known vectors v1and v2, get vectors d1and d2by rotation and scaling. (c) Get vector das d2−d1. Reject either d1 or d2on dto compute p. the desired vector pis the rejection of d1or d2on d(see Fig. 2c). In VGA, the above steps are summarised in the following equation p=(d1∧d)d−1=−(d2∧d)d−1=(d1∧d2)d−1(16) It should be noted that the proposed solution is remarkably simple and does not involve the use of any type of coordinates. The result is obtained by simple geometric operations, such as rotation, scaling, and rejection, applied to the inherent primitives of VGA, such as vectors in this specific case. The proposed method has advantages over the existing ones. It avoids some limitations as in Tienstra’s method where no solution can be found if the points A,B, and C are collinear. Furthermore, it is feasible to obtain an indicator of how close Pis to the forbidden circle (defined by points A,B, and C) by means of the length of vector d.If the point Pis on this circle, then it can be easily checked that d=d=0. Consequently, small values of dsuggest that we should relocate the station to another site to ensure a reduced source of error (see Sect. 6for a detailed error analysis). 4.2 Conformal GA method The methodology employing VGA, as outlined in Sect. 4.1, primarily utilised the GA G2. However, given that the resection problem primarily deals with circles, utilising its conformal extension, specifically the compass ruler algebra (CRA), appears to be more suitable (see Hildenbrand 2018). When dealing with geometric problems, it is highly advantageous to have a tool that can express graphical methods algebraically. On the basis of the postulates presented in Sect. 2.2,two traditional and well-known graphical methods are proposed to solve the resection problem using CRA: Cassini and Collins. By leveraging CRA, both methods can receive clear algebraic interpretations, as explained in the following sections. For a more in-depth example of CRA using the Clifford library in Python, see the GitHub repository. 4.2.1 Cassini construction The Cassini method provides a solution to the resection problem by leveraging the inscribed angle theorem. The solution is obtained by determining the intersection of two circles: one passing through points A,B, and P, and the other through points B,C, and Pas shown in Fig. 3. To determine the centres of the circles, two lines must be intersected. The stepby-step graphical approach underlying the Cassini method can be elucidated, along with the equivalent steps, using the CRA algebra (see Fig. 4). 123
47 Page 8 of 21 J. Ventura et al. Fig. 3 2D graphical resection procedure using Cassini method 1. CRA Mapping: The problem starts by mapping three known Euclidean points, A,B, and C, to the CRA domain. For example, the point Awith coordinates (a1,a2)is mapped as A=(a1,a2)−→ a=a1σ1+a2σ2 −→ a=σ0+ a+1 2 a2σ∞(17) 2. Auxiliary Lines: To obtain the position of the centre of the circle defined by A,B, and P, the line AB joining A and Bmust first be constructed along with the line mAB as the perpendicular bisector of AB. The line AO is then constructed by rotating AB by (π/2−α) relative to A clockwise. This process is repeated for BC and mBC but with CO being rotated some (π/2−β)counterclockwise relative to C. LAB =a∧b∧σ∞,LCB =c∧b∧σ∞ MAB =(a−b)MBC =(b−c)(18) For the rotated lines AO1and CO 2, the following rotors and translators must be first defined Rα=e−1 2(α−π/2)σ12 ,Rβ=e−1 2(π/2−β)σ12 TA=1− a 2 σ∞,TC=1− c 2 σ∞ DA=TARα TA,DC=TCRβ TC(19) , and thus, LAO1=DALAB DA,LCO 2=DCLBC DC(20) Note that rotations are always relative to the origin. To rotate around an arbitrary point, it must first be translated to the origin. The rotation is then applied, followed by another translation that returns the point back to its original position. 3. Circle building: The centres of the circles can now be found as the intersection of the lines AB and AO1and BC and CO 2 O1=o1∧σ∞=LAB ∨LAO =(L AB ∧L AO) O2=o2∧σ∞=LBC ∨LCO =(L BC ∧L CO)(21) Fig. 4 Graphical solution of Cassini method step by step. The point Pis found by intersecting the two circles c1 and c2where all the points see AB with angle αand BC with angle β, respectively (a) Lines defined by AB and BC (b) Bisector of AB and BC (c) Centres O1and O2(d) The sought point P 123
A novel geometric method based on conformal geometric algebra applied to the resection problem… Page 9 of 21 47 Fig. 5 CRA version of Cassini’s method step by step The result is a flat point, representing the wedge of the sought point and the point at infinity.1The extraction of the point of interest is straightforward by factoring out the point at infinity o1=(σ0·O1−σ0) o2=(σ0·O2−σ0)(22) The radius of the circles can be computed as r1=−2(o1·a), r2=−2(o2·c)(23) , and the circles themselves are determined as c1=o1−1 2r2 1σ∞c2=o2−1 2r2 2σ∞(24) 4. Intersection of circles: The desired result can be obtained from the intersection of the two circles c1and c2as P=c1∧c2(25) Two intersecting circles in CGA yield a couple of points (1D sphere), also known as pair-point P. Finally, the 1A line can be considered as a circle with infinite radius, so the intersection of two lines results in two points, one at infinity. Euclidean point Pis recovered by the classic formula (see Hildenbrand 2018) P±=±√P·P+P σ∞·P(26) One of the points P±is exactly the point B, and the other is the sought point P(see Fig. 3). Figure 5shows a concise and condensed summary of the key steps involved and discussed above. It provides valuable visual depictions that enhance the geometric intuition underlying the method, offering an algebraic interpretation of the graphical approach. It also helps to reinforce the strong connection between the graphical Cassini method and its algebraic translation using CRA. The GitHub repository shows several examples derived from the code developed by the authors. The computational procedures in CGA are found to be straightforward. The methodology involves the manipulation and combination of geometric objects, thus justifying the occasional reference to GA as an object-oriented approach. 4.2.2 Collins construction The graphical method of Collins provides a solution to the resection problem using the intersection of the line passing through the point Band the so-called Collins auxiliary point 123
47 Page 16 of 21 J. Ventura et al. Table 5 Performance comparison of resection algorithms using geometric algebra and state-of-the-art algorithms (see Pierlot and Van Droogenbroeck 2014) Method Mean (μs) Rank # VGA 18997.8 1 Total #1 21919.1 2 Total #2 27340.0 3 CollinsCGA 196122.2 4 CassiniCGA 220135.8 5 Our findings reveal that our VGA-based algorithm outperforms state-of-the-art methods (see Table 5), executing approximately 13.3% and 30% faster than the previously best-known algorithms by Pierlot with Total #1 and Total #2. The primary advantage of our methodology and implementation is the utilisation of GA-FuL’s comprehensive code generation capabilities. These capabilities range from generating code for individual multivector operations to creating full software libraries with proper software architecture and nested folder/file structure, enabling efficient and optimised geometric algebra computations. 6 Uncertainty analysis This section investigates the impact of measurement uncertainties on the efficacy of the proposed GA methods. Given the intrinsic presence of noise in practical measurements, it is crucial to assess the sensitivity and resilience of the method to such perturbations (Font-Llagunes and Batlle 2009a;Pierlot and Van Droogenbroeck 2014). Previous discussions have highlighted that the resection problem faces an intractable challenge within the forbidden circle defined by points A,B, and C. This limitation is inherent in the nature of the problem rather than a reflection of the methods used to solve it. Novel methodological approaches should be considered with caution as they may inadvertently introduce undesirable scenarios and singularities not covered by the original formulation of the problem. In this regard, the proposed methods handle degenerate or limiting cases in a convenient and elegant graphical way, as shown in Fig. 11. This study presents three approaches to address the problem of 2D and 3D resection. In order to improve the ability to evaluate the accuracy of the algorithms and determine the station position error, a new metric has been developed. Consequently, we propose three formulations that are essentially different variations of the same underlying approach, each defining the metric Dbased on the square of a distance to denote the proximity of the station Pto the forbidden region. •VGA method: The vector dis used as the difference between the vectors d1and d2(refer to Eq. (15)). When the station Pis located on the forbidden circle, the vectors d1and d2converge, resulting in a null vector d. Therefore, minimal values of Dindicate close proximity to the forbidden circle D=d2(41) •Cassini method: The determination of Dusing the CGA Cassini method requires the calculation of the squared distance between the centres of the circles c1and c2(as shown in Fig. 4), marked by O1and O2, respectively. In situations where the station is positioned on the forbidden circle, the two circles converge, causing their centres to overlap and the intercentre distance to become 0. Thus, according to Eq. (A3), Dis defined as D=−2(o1·o2)(42) •Collins method: Similarly to the previous method, Dis determined as the squared distance between two specific points. The auxiliary Collins point Eand the point Bare used for this purpose. They coincide when the station P is situated on the forbidden circle. The formula used to compute this distance is as follows D=−2(e·b)(43) To validate the sensitivity of the proposed algorithms, a series of simulations have been proposed. The simulation framework is designed within a square area measuring 4 by 4 square metres, incorporating two unique configurations for three known points. The first configuration is an equilateral triangle, with the points located at the positions A=(0,1), B=(−0.866,−0.5), and C=(0.866,−0.5). The points in the second configuration are linearly arranged at A= (−0.866,0),B=(0,0), and C=(0.866,0). A spacing of 2cm is used in each direction across the grid. At each point on the grid, the angles αand βseen from Pare calculated. Gaussian noise is introduced into these angles, characterised by a zero mean and two distinct standard deviations (σ= 0.01 degrees and σ= 0.1 degrees). The algorithms use these modified angles as input to determine the estimated position of the unknown point. The discrepancy in position (d)is quantified by the Euclidean distance between the exact and estimated location of the point P. The study performs 1000 iterations for each position to determine the standard deviation of the position error. The resulting standard deviations are shown in Fig.14. The results for the equilateral triangle configuration, with standard deviations of σ= 0.01 degrees and σ= 0.1 degrees, are presented in the first and second columns, respectively. Similarly, the 123
A novel geometric method based on conformal geometric algebra applied to the resection problem… Page 17 of 21 47 Fig. 14 Error analysis for two point configurations using a square window of 2 ×2 m.. Config #1 is an equilateral triangle, while in config #2 the points are collinear. Both configurations have been tested with standard deviation σ=0.01 and σ=0.1. The three GA methods (VGA, CollinsCGA, and CassiniCGA) are compared based on their position error (first row) and metric 1/D(second row), respectively 123
47 Page 18 of 21 J. Ventura et al. results for the second configuration (collinear points) are detailed in the third and fourth columns, respectively, using the same standard deviations. The figure shows the standard deviation of the position error and the mean error measure 1/Din the first and second rows, in that order. It is important to note that the scales used in the graphic representation are not linear. To enhance the visual clarity of these images and emphasise the correlation between the position errors and our new error metric, we applied histogram equalisation to the images. 6.1 Discussion of the results The simulations performed are in agreement with those reported in the literature and support the case studies of each of the three methods presented. Figure 14 shows the different point configurations, where the forbidden circle is clearly identifiable due to the increase in standard deviation of the position error as it is approached. Minimal errors are observed inside the circle, while errors increase with distance outside the circle. It is remarkable that the three methods provide almost identical results for the error position. The 1/Dmetric plots maintain a shape similar to the position error plots near the forbidden circle. High red values are observed as Dbecomes smaller, while minimum values are observed in the lines AB and BC in Cassini and AC in Collins due to an increase of D(see Figs. 9and 10). Configuration 2 transforms the forbidden circle into a straight line. In this region, the error distribution is high and increases progressively with distance from the line. Because this configuration represents a degenerate situation, the case where Dbecomes zero does not occur. The different methods used to define Dlead to overlapping points at infinity, causing Dto take high values near the critical line and resulting in a minimal value for 1/D. The metric Dwas chosen for configuration 2 due to this property. All methods produce identical position error plots, indicating consistency and conformity with the results obtained by most algorithms to solve the resection problem. This confirms that the sensitivity to calculate the position of the point P, even with noisy measured angles, is independent of the method used and unique, as discussed in Pierlot and Van Droogenbroeck (2014) and Font-Llagunes and Batlle (2009b). The metric 1/Dcan be used as an indicator of proximity to the forbidden circle. When dealing with aligned beacons, the value of Dshould be used directly. In other cases, the similarity between this metric and the position error suggests that the former can approximate the latter if a function of the other problem parameters is applied. 7 Conclusions This article presents a novel approach to solving the resection problem in two and three dimensions using conformal geometric algebra (CGA). The CGA framework allowed representing points, lines, planes, circles, and spheres in a unified mathematical structure and offered a more intuitive understanding and efficient solution to the resection problem compared to existing algebraic techniques. The proposed method leveraged the ability of CGA to transition between different reference frames without requiring coordinate transformations. This eliminated the need for multiple calculation steps and complex algebraic manipulations that are characteristic of traditional algebraic solutions. Through extensive numerical simulations, we have demonstrated the validity and efficacy of our GA-based approach, achieving accuracy comparable to that of established algebraic techniques, while significantly improving computational efficiency and providing valuable geometric insights. Our findings suggest that the geometric algebra framework has strong potential for solving resection-type problems not only in surveying and geodesy but also in computer graphics, robotics, computer vision, and navigation. By exploiting geometric relationships between entities, CGA paves the way for more intuitive solutions that unify computations involving different geometric primitives. Future research can build upon the ideas presented here to address more complex variants of the resection problem that involve additional constraints. The CGA method can also be extended to address intersection problems and other spatial geometric computations across diverse disciplines. By harnessing the power of geometric algebra and the versatility of conformal geometric methods, this work opens up new possibilities for advancing geometric research and computational techniques. Appendix A Some calculations over CGA CGA has been briefly introduced in Sect.2.2.Inthis Appendix, some fundamental calculations are presented to illustrate the computational efficiency inherent in this algebra. For further information and a detailed understanding, the reader is referred to Hrdina et al. (2021), Dorst et al. (2010), Hestenes and Sobczyk (2012), Hildenbrand (2018). A.1 Conformal inclusion We’ll show how the translations and rotations are connected and related to the insertion of the Euclidean space. If we choose vector σ0as the origin of the coordinate system and use the element 123
A novel geometric method based on conformal geometric algebra applied to the resection problem… Page 19 of 21 47 T=e1 2(xσ1+yσ2+zσ3)∧σ∞=1 +1 2(xσ1+yσ2+zσ3)∧σ∞(A1) , then it can act as a translation, allowing us to define an inclusion as a translation of the origin σ0in the direction of the vector x=(xσ1+yσ2+zσ3)by the sandwich product as ι( x)= Tσ0T=1−1 2 x∧σ∞σ01+1 2 x∧σ∞ =σ0−1 2 xσ∞σ01+1 2 xσ∞ =σ0+σ0 1 2 xσ∞−1 2 xσ∞σ01+1 2 xσ∞ =σ0+σ0 1 2 xσ∞−1 2 xσ∞σ0−1 2 xσ∞σ0 1 2 xσ∞ =σ0−1 2 x(σ0σ∞+σ∞σ0) − x21 4 σ∞σ0σ∞=σ0+ x+1 2 x2σ∞ So, we see that we can identify points in R3with vectors in CGA through the inclusion ι: x=(x,y,z)→ x=σ0+ x+1 2 x2σ∞(A2) where the element Tacts as a translation. Using a similar reasoning, the element R=en1σ2σ3+n2σ1σ3+n3σ1σ2acts as a rotation due to identification Im H=σ2σ3,σ1σ3,σ1σ2 and the following computation Rι( x) R=Rx R=Rσ0+ x+1 2 x2σ∞ R =Rσ0 R+R x R+1 2 x2Rσ∞ R =σ0+R x R+1 2 x2σ∞ A.2 Distance between points One of the key distinctions between CGA and VGA lies in the interpretation of the inner product. In CGA, the inner product represents the distance between points, which is fundamentally a quadratic concept in VGA. By incorporating two additional dimensions, CGA linearises this quadratic function. The equation below elucidates this: a·b=σ0+ a+1 2 a2σ∞·σ0+ b+1 2 b2σ∞ = a· b−1 2 a2−1 2 b2=−1 2( a− b)2(A3) This calculation explicitly demonstrates how CGA linearises the quadratic nature of the distance representation found in VGA, offering nuanced insights into geometric relations and interactions between points. Furthermore, this characteristic is utilised to depict the midline Mas illustrated in expression (18). A midline is established by two control points, aand b, and encompasses points that maintain equal distances to these control points. This means that (x−a)2=(x−b)2. In the context of CGA, this condition can be articulated as follows: x·a=x·b⇒0=x·a−x·b=x·(a−b), i.e. x∈MAB ↔x·(a−b)=0, so (a−b)is the IPNS representation of the midline defined by control points aand b. A.3 IPNS vs OPNS representation of objects In GA, we have several types of product, so it is possible to represent objects in different ways. With the help of the inner product, we can define the object Cas follows. x∈C⇔x·C=0(A4) In CGA, for example, the object C=n1σ1+n2σ2+n3σ3+ dσ∞can be tested and verify that it represents a plane (σ0+ x+1 2 x2σ∞)·( n+dσ∞)= x· n−d=0(A5) where n=n1σ1+n2σ2+n3σ3. Equation (A5) describes a plane with normal vector nand distance from the origin d. On the other hand, the wedge product defines an object Das follows x∈D⇔x∧D=0. Again, in CGA, for example, the object D=a∧b∧σ∞=σ0+ a+1 2 a2σ∞ ∧σ0+ b+1 2 b2σ∞∧σ∞ =(σ0+ a)∧(σ0+ b)∧σ∞ =σ0∧ b∧σ∞+ a∧σ0∧σ∞+ a∧ b∧σ∞ defines the line goes through the points aand b: x∧D=σ0+ x+1 2 x2σ∞ ∧(σ0∧ b∧σ∞+ a∧σ0∧σ∞+ a∧ b∧σ∞) 123
47 Page 20 of 21 J. Ventura et al. =σ0∧ a∧ b∧σ∞+ x ∧(σ0∧ b∧σ∞+ a∧σ0∧σ∞+ a∧ b∧σ∞) =(σ0∧( a∧ b+ x∧( a− b)) + x∧ a∧ b)∧σ∞ =0 So a∧ b+ x∧( a− b)=0 and x∧ a∧ b=0,(A6) which represented the line based on the points aand b. Funding Funding for open access publishing: Universidad de Almería/CBUA. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/. References Awange JL, Grafarend EW (2005) Solving algebraic computational problems in geodesy and geoinformatics. Algebraic computational problems in geodesy and geoinformatics Awange JL (2002) Groebner basis solution of planar resection. Surv Rev 36(283):528–543 Awange JL, Grafarend EW (2002) Sylvester resultant solution of planar ranging problem. Allgemeine Vermessungs-Nachrichten 108(4):143–146 Awange J, Grafarend EW (2003a) Groebner-basis solution of the threedimensional resection problem (p4p). J Geodesy 77:327–337 Awange JL, Grafarend EW (2003b) Multipolynomial resultant solution of the three-dimensional resection problem (p4p). Bollettino di Geodesia e Science affini 62(2):79–102 Awange JL, Grafarend EW, Paláncz B, Zaletnyik P (2010) Algebraic geodesy and geoinformatics. Springer, Berlin Bähr H (1991) Einfach überbestimmtes ebenes einschneiden, differentialgeometrisch analysiert. Zeitschrift für Vermessungswesen 116(1991):545–552 Bil P (1992) Sectie en projectie. Nederlands Geodetisch Tijdschrift Geodesia Bock W (1956) Mathematische und geschichtliche betrachtungen zum einschneiden. Ph.D. thesis, Institute f. Geodäsie u. Photogrammetrie d. Technischen Hochschule Cameron J, Lasenby J (2008) Oriented conformal geometric algebra. Adv Appl Clifford Algebras 18:523–538 Dekov D (2012) A numerical method for solving the horizontal resection problem in surveying. J Geodetic Sci 2(1):65–67 Dorst L, Fontijne D, Mann S (2010) Geometric algebra for computer science: an object-oriented approach to geometry. Elsevier, Burlington Eid AH Geometric Algebra Fulcrum Library (GA-FuL). https://github. com/ga-explorer/GeometricAlgebraFulcrumLib. Accessed 2024 Fischler MA, Bolles RC (1981) Random sample consensus: a paradigm for model fitting with applications to image analysis and automated cartography. Commun ACM 24(6):381–395 Font-Llagunes JM, Batlle JA (2009a) Consistent triangulation for mobile robot localization using discontinuous angular measurements. Robot Auton Syst 57(9):931–942 Font-Llagunes JM, Batlle JA (2009b) New method that solves the threepoint resection problem using straight lines intersection. J Surv Eng 135(2):39–45 Grunert JA (1841) Das pothenot’sche problem, in erweiterter gestalt nebst bemerkungen über seine anwendung in der geodäsie". Archiv der Mathematik und Physik 1:238–248 Hadfield H, Wieser E, Arsenovic A, Kern R (2021) The Pygae Team: Pygae/clifford. https://doi.org/10.5281/zenodo.1453978 Hestenes D, Sobczyk G (2012) Clifford algebra to geometric calculus: a unified language for mathematics and physics, vol 5. Springer, Dordrecht Hildenbrand D (2018) Introduction to geometric algebra computing, 1st edn. Chapman and Hall/CRC, Boca Raton Hitzer E, Lavor C, Hildenbrand D (2022) Current survey of clifford geometric algebra applications. Math Methods Appl Sci Hrdina J, Návrat A (2017) Binocular computer vision based on conformal geometric algebra. Adv Appl Clifford Algebras 27:1945–1959 Hrdina J, Návrat A, Vašík P, Matoušek R (2017) CGA-based robotic snake control. Adv Appl Clifford Algebras 27:621–632 Hrdina J, Návrat A, Vašík P, Dorst L (2021) Projective geometric algebra as a subalgebra of conformal geometric algebra. Adv Appl Clifford Algebras 31:1–14 Masselli A, Zell A (2014) A new geometric approach for faster solving the perspective-three-point problem. In: 2014 22nd international conference on pattern recognition, pp 2119–2124. IEEE Mazaheri M, Habib A (2015) Quaternion-based solutions for the single photo resection problem. Photogramm Eng Remote Sens 81(3):209–217 McCaw GT (1918) Resection in survey. Geograph J 52(2):105–123 Montoya FG, Baños R, Alcayde A, Arrabal-Campos FM (2019) Analysis of power flow under non-sinusoidal conditions in the presence of harmonics and interharmonics using geometric algebra. Int J Electr Power Energy Syst 111:486–492 Montoya FG, Baños R, Alcayde A, Arrabal-Campos FM, Roldán-Pérez J (2021) Vector geometric algebra in power systems: An updated formulation of apparent power under non-sinusoidal conditions. Mathematics 9(11):1295 Paláncz B, Awange JL, Zaletnyik P, Lewis RH (2010) Linear homotopy solution of nonlinear systems of equations in geodesy. J Geodesy 84:79–95 Pereira FI, Luft JA, Ilha G, Susin A (2018) A novel resectionintersection algorithm with fast triangulation applied to monocular visual odometry. IEEE Trans Intell Transp Syst 19(11):3584–3593 Pierlot V, Van Droogenbroeck M (2014) A new three object triangulation algorithm for mobile robot positioning. IEEE Trans Robot 30(3):566–577 Selig JM (2005) Geometric fundamentals of robotics, 2nd edn. Monographs in computer science. Springer, New York. https://doi.org/ 10.1007/b138859 Smith J (2023a) Solving the Snellius-Pothenot resection (surveying) problem via geometric algebra. https://www.youtube.com/watch? v=h863AAQ3lF8. Accessed 9 April 2023 Smith J (2023b) Via geometric algebra: a solution to the SnelliusPothenot resection (surveying) problem. arXiv:2305.0079 Sturmfels B (2002) Solving systems of polynomial equations, vol 97. American Mathematical Society, Berkeley Wareham R, Cameron J, Lasenby J (2004) Applications of conformal geometric algebra in computer vision and graphics. In: Inter123
A novel geometric method based on conformal geometric algebra applied to the resection problem… Page 21 of 21 47 national workshop on mathematics mechanization, pp 329–349. Springer, Berlin Wreede LC (2007) Willebrord Snellius (1580–1626): a Humanist Reshaping the Mathematical Sciences. Utrecht University, Utrecht Zaplana I, Hadfield H, Lasenby J (2022) Closed-form solutions for the inverse kinematics of serial robots using conformal geometric algebra. Mech Mach Theory 173:104835 123