Orthogonal distance from an ellipsoid

Distancia ortogonal a partir de um elipsóide

Abstracts

Finding the orthogonal (shortest) distance to an ellipsoid corresponds to the ellipsoidal height in Geodesy. Despite that the commonly used Earth reference systems, like WGS-84, are based on rotational ellipsoids, there have also been over the course of the years permanent scientific investigations undertaken into different aspects of the triaxial ellipsoid. Geodetic research has traditionally been motivated by the need to approximate closer and closer the physical reality. Several investigations have shown that the earth is approximated better by a triaxial ellipsoid rather than a rotational one Burša and Šima (1980). The problem of finding the shortest distance is encountered frequently in the Cartesian- Geodetic coordinate transformation, optimization problem, fitting ellipsoid, image processing, face recognition, computer games, and so on. We have chosen a triaxial ellipsoid for the reason that it possesess a general surface. Thus, the minimum distance from rotational ellipsoid and sphere is found with the same algorithm. This study deals with the computation of the shortest distance from a point to a triaxial ellipsoid.

Orthogonal (Shortest) Distance; Triaxial Ellipsoid Coordinate Transformation; Fitting Ellipsoid

Encontrar a distância orthogonal a um elipsóide corresponde a altura elipsoidal em Geodésia. Apesar de os sistemas de referência da Terra mais comumente usados, como WGS-84, são baseados em elipsóides rotacionais, tem tido por anos, investigações científicas permanentes feitas em diferentes aspectos do elipsóide triaxial. A pesquisa geodésica tem sido tradicionalmente motivada pela necessidade de uma aproximação cada vez mais próxima da realidade física. Diversas investigações, tem mostrado que a Terra se aproxima mais de um elipsóide traxial ao invés de rotacional, Busa e Sima (1980). O problema de se encontrar a distância mais curta é encontrado frequentemente em transformações de coordenadas geodésicas cartesianas, problema de otimização, adequação de elipsóide, processamento de imagens, reconhecimento de faces, jogos de computador, etc. Nós escolhemos o elipsóide triaxial pela razão de que ele tem uma superfície geral. Assim, a distência mínima a um elipsóide rotacional e a esfera é encontrada com o mesmo algorítimo. Este trabalho realiza o cálculo da distância mais curta a partir de um ponto de um elipsóide triaxial.

Distância Ortogonal; Elipsóide Triaxial; Transformação de Coordenadas

1. INTRODUCTION

Although ellipsoid (general or triaxial) equation Eq. (1) is quite simple and smooth but geodetic computations are quite difficult on the ellipsoid. The main reason for this difficulty is the lack of symmetry. Triaxial ellipsoid is generally not used in geodetic applications. Rotational ellipsoid (ellipsoid revolution, biaxial ellipsoid, spheroid) is frequently used in geodetic applications.

Today increasing GPS and satellite measurement precision will allow us to determine more realistic earth ellipsoid. Geodetic research has traditionally been motivated by the need to continually improve approximations of physical reality.

Geodetic research has traditionally been motivated by the need to approximate closer and closer the physical reality. Several investigations have shown that the earth is approximated better by a triaxial ellipsoid rather than a rotational one.

Furthermore, non-spherical celestial bodies such as planets, physical satellites, asteroids and comets can be modeled by a triaxial ellipsoid. Also, the present day accuracy requirements and the modern computational capabilities push toward the study on the triaxial ellipsoid as a geometrical and a physical model in geodesy and related interdisciplinary sciences Panou et al(2013).

First, the basic definition of ellipsoid starts with giving mathematical equations to explain the concepts. To show how computations of the shortest distance to an ellipsoid, are carried out, we solve this problem separately: standard ellipsoid and the shifted-oriented ellipsoid. The efficacy of the new algorithms is demonstrated through simulations.

When we look at the literature on this subject we see the various studies: Eberly (2008), Feltens (2009) and Ligas (2012). For the solution Eberly (2008) gives a method that is based on sixth degree polynominal. He has benefited from the largest root of 6th degree polynomial. Feltens (2009) gives a vector-based iteration process for finding the point on the ellipsoid. Ligas (2012) claims his method turns out to be more accurate, faster and applicable than Feltens method. The presented paper tries to give the shortest distance calculation not only for the ellipsoid in standard position but also the shifted-oriented ellipsoid. We could not find enough studies with numerical examples on this subject in the literature, especially for shifted-oriented ellipsoid.

Triaxial ellipsoid formulas are quite useful, because obtaining the rotational ellipsoid formula from triaxial ellipsoid formula is easy. For this, equatorial semi-axis are accepted equal to each other (a_x=a_y =a) which is sufficient on triaxial ellipsoid formula. Similarly to obtain sphere formula from rotational ellipsoid formula it is sufficient to take as (a = b= R) Bektas (2009).

1.1 Ellipsoid

An ellipsoid is a closed quadric surface that is analogue of an ellipse (see Figure 1).

Figure 1
- The shortest distance on triaxial ellipsoid (standard position)

Ellipsoid has three different axes (ax>ay>b). Mathematical literature often uses "ellipsoid" in place of "Triaxial ellipsoid or general ellipsoid". Scientific literatura (particularly geodesy) often uses "ellipsoid" in place of "biaxial ellipsoid, rotational ellipsoid or ellipsoid revolution". Older literature uses 'spheroid' in place of rotational ellipsoid. The standard equation of an ellipsoid centered at the origin of a cartesian coordinate system and aligned with the axes is given below. (http://en.wikipedia.org/wiki/Ellipsoid).

2. FINDING THE POINT ON THE ELLIPSOID

It can be proved that the shortest distance is along the surface normal. The first step is to find the projection of an external point denoted as P_G(x_G,y_G,,z_G ) as shown in Figure 1 onto this ellipsoid along the normal to this surface i.e. point P_E (x_E, y_E,,z_E) Feltens (2009), Ligas (2012). P_G (x_G, y_G,,z_G ) is a point on the earth surface. This section is handled under two headings: first, the shortest distance from the standard ellipsoid, second, the shortest distance from the shifted-oriented ellipsoid.

2.1. The Shortest Distance from the Standard Ellipsoid

In this section, we assume that the ellipsoid is in standard position, in other words its axis is aligned and centered at the origin (Figure 1). We will discuss the other ellipsoid under the next title. The following definitions will be used.

ax = equatorial semimajor axis of the ellipse ay = equatorial semiminor axis of the ellipse b = polar semi-minor axis of the ellipse λ = geodetic longitude ψ = geodetic latitude h = PG PE : ellipsoid height : the shortest distance

Feltens (2009) gives a vector-based iteration process for finding the point on the ellipsoid. Ligas (2012) claims his methods turn out to be more accurate, faster and applicable. Ligas' Method is based on solving nonlinear system of equation. Ligas' Method is as follows:

In accordance with Figure 1 a collinearity condition can be written between P_E and P_G

The first stage begins with constructing two collinear vectors: a vector normal (n) to the ellipsoid (obtained from the gradient operator of a triaxial ellipsoid) in the point P_E that may be expressed as (seen in Figure 1):

From the above Eq. (8) Ligas wrote the below given three equations:

But only two of them are necessary to set a single variant of a system of nonlinear equations to be solved. In addition, the coordinates of P_E must satisfy the equation of the triaxial ellipsoid i.e.

From Eqs(9.a),(9. b) and (9.c) three pairs of equations are obtained which together with fixed one Eq.(9.d) produce three variants of the system of equations to be solved in order to obtain the solution for P_E. As a result, it is established three nonlinear sytems of equation:

Case.1:Eqs.(9.a),(9. b) and (9.d) Case.2:Eqs.(9.a),(9. c) and (9.d) Case.3:Eqs.(9.b),(9. c) and (9.d)

Ligas (2012) has been doing research on three different nonlinear system of equation in terms of run-time and the number of iterations. But in our opinion it is pointless. It is not possible to say one of the systems will yield better results than the others. The distinctions between the results of different systems are actually meaningless. The distinctions completely arise from rounding errors. Because only two of the three equations ((9.a), (9. b) and (9.c)) are independent. Hence we can choose any two of these three equations. For example, Case.1 is chosen and these three equations are linearized by Taylor series expansion and the system of equations is solved in order to obtain the solution for X_E = [x_E , y_E , z_E ].

The initial guesses for the point on the ellipsoid P_E were chosen the same as in Feltens (2009), namely:

For the solution of nonlinear system (9.a,9.b,9.d) to obtain the value for X_E = [x_E , y_E , z_E] is necessary the calculus of the Jacobi's Matrix , the vector f, calculated in equations (11.a,11.b,11.c).

The entries may be written as follows

If δ_Ei is less than threshold, the iteration is stopped. After the first step is accomplished, finding P_E on the ellipsoid, and having coordinates of P_E , the shortest distance P_G P_E,= h may easily be computed as:

The following link can be used for the shortest distance and projection coordinates on triaxial ellipsoid Bektas (2014). http://www.mathworks.com/matlabcentral/fileexchange/46261-the-shortest-distance-from-a-point-to-ellipsoid.

2.2. The shortest distance from the shifted-oriented ellipsoid

In general, for shifted-oriented ellipsoid as in (Figure 2), the data point P_G(X_G,Y_G,,Z_G ) can be rotated and translated to axis-aligned ellipsoid centered at the origin and the distances can be calculated in that system.

Figure 2
- Shifted - oriented ellipsoid

For this conversion we utilize performed as follows by making use of ellipsoid's center coordinates (Xo,Yo,Zo) and the rotation angles (ε, ψ, ω) , in accordance with Figure 2

The shortest distance ( h ) calculation is made of the new converted coordinates P_G(x_G, y_G,z_G) in standard position with the above procedure and the coordinates of P_E (x_E,y_E,,z_E) in standard position are found. To find the true coordinates of P_E(X_E, Y_E, Z_E) we need to make a transformation as below:

3. Obtaining Ellipsoid Parameters from Conic Equation

Generally an ellipsoid is defined with 9 parameters. These parameters are; 3 coordinates of centre (X_o,Y_o,Z_o), 3 semi-axes (a_x,a_y,b) and 3 rotational angles (ε, ψ, ω) which represent rotations around x-,y- and z- axes respectively as shown in Figure 2. These angles control the orientation of the ellipsoid.

If we know the conic equation, ellipsoid parameters can be calculated. Here is a conic equation

Necessary conditions for this problem to have a unique solution are that n > = 9, (n: denotes the number of data point cartesian coordinates(x,y,z)) and the data points lie in general position (e.g., do not all lie is an elliptic plane).Throughout this paper, we assume that these conditions are satisfied. v = [ A B C D E F G H I ]^Tunknown conic parameters.

It is solved easily in the Least Square (LS) sense by MATLAB as below:

If n=9 the solution is named exact solution. All data points are satisfy Eq. (21) .

If n>9 all data points do not satisfy Eq.(21). The solution named is algebraic ellipsoid fitting.

In theory, the conditions that ensure a quadratic surface to be an ellipsoid have been well investigated and explicitly stated in analytic geometry textbooks. Since an ellipsoid can be degenerated into other kinds of elliptic quadrics, such as an elliptic paraboloid. Therefore a proper constraint must be added by Li and Griffiths (2004) gives the following definitions:

However 4j-i2 > 0 is just a sufficient condition to guarantee that an equation of second degree in three variables represent an ellipsoid, but it is not neccesary. In this paper, we assume that these conditions are satisfied.

After conical equation Eq. (21) is solved in the LS sense by Eq. (22), this section we determine the center, semi-axis and rotation angles of the ellipsoid using an algorithm from Yury Petrov's Ellipsoid Fit Method Petrov (2009). MATLAB script the following link http://www.mathworks.com/ matlabcentral/fileexchange/24693-ellipsoid-fit.

The solution of the above equation system which is established with conic coefficients gives us the coordinates of ellipsoid's center (X_o,Y_o,Z_o).

In order to find of semi-axis (a_x,a_y,b) and rotation angles (ε, ψ, ω) of the ellipsoid, we use the following M and T matrix,

Eigenvalues and eigenvectors of this matrix S can be easily calculated whit the following MATLAB command,

evals : Eigenvalues of Let's assume λ1 , λ2 and λ3 are the eigenvalues of the matrix S, in descending order.

Semi-axes of ellipsoid (ax,ay,b) are obtained from the eigenvalues of S as below;

3.1 Numeric Example-1 (for standard ellipsoid)

In order to demonstrate the validity of the shortest distance algorithms presented above, a numerical example is given. The algorithm was implemented in MATLAB. The numerical computations in the triaxial case were carried out using Earth's geometrical parameters a_x= 6378388.0000m, a_y = 6378318.0000m and b= 6356911.9461m.

The following link can be used for shortest distance from a point to triaxial ellipsoid http://www.mathworks.com/matlabcentral/fileexchange/46261-the-shortest-distance-from-a-point-to-ellipsoid.

Given P_G point Cartesian coordinates x_G = 3909863.9271m, y_G =3909778.1230 m, z_G = 3170932.5016m. Finding the shortest distance to triaxial ellipsoid.

For this, we must firstly find P_E point Cartesian coordinates. The initial guesses for the P_E point on the ellipsoid are determined according to Eq.(10)

Finding P_E point on the ellipsoid from P_G point iteratively from Eqs.(12-16) ;

3.2 Numeric Example-2(for shifted-oriented ellipsoid)

12 triplets (x,y,z) cartesian coordinates are given below.

Here are data points coordinates:

a) Find the conic equation and the ellipsoid parameter; b) Given P_G point cartesian coordinates P_G (X_G = 7 , Y_G =22 , Z_G = 31), find P_E point cartesian coordinates and the shortest distance from the ellipsoid.

Solution

a) The conical coefficents in the LS sense from Eq.(22) is,

v = [-0.0006 -0.0008 -0.0010 0.0005 -0.0005 0.0003 0.0092 0.0050 0.0278]

We determine the center, semi-axis and rotation angles of the ellipsoid using an algorithm from Yury Petrov's ellipsoid fit method.

The result of ellipsoid parameter is shown in Table 1:

Table 1
- The result of ellipsoid parameter

	Center of Coordinates			Semi-axes							Rotate angles
												(degree)
	Xo	Yo	Zo	ax		ay			b		ε		ψ	ω
	10.3837	20.9653	29.0070	7.4676		3.1643		2.0147		47.98			18.68	28.21

b) P_G (X_G = 7 , Y_G =22 , Z_G = 31)

The projection coordinates P_E (x_E, y_E,,z_E)

The following link can be used for projection coordinates on triaxial ellipsoid http://www.mathworks.com/matlabcentral/fileexchange/46261-the-shortest-distance-from-a-point-to-ellipsoid

Here, minus sign indicates that point P_G remains in the ellipsoid

4. CONCLUSION

In this study, we handled the computation of the shortest distance from a point to an ellipsoid. We have shown this in two separate cases (standard ellipsoid and shifted-oriented ellipsoid). The problem of finding the shortest distance is encountered frequently in the Cartesian- Geodetic coordinate transformation, optimization problem, fitting ellipsoid, image processing, face recognition, computer games, and so on. In conclusion, the presented method may be considered as fast, accurate and reliable and may be successfully used in other areas. The presented algorithm can be applied easily for rotational ellipsoid and sphere and, also for other surface such as paraboloid, hyperboloid.

REFERENCES

• BEKTAS, S. "Practical Geodesy", Ondokuz Mayis University Press,ISBN 978-975-7636-65-6. . (2009).
• BEKTAS, S. The shortest distance from a point to Ellipsoid. Access on: 23/07/2014. http://www.mathworks.com/matlabcentral/fileexchange/46261-the-shortest-distance-from-a-point-to-ellipsoid. (2014).
  » http://www.mathworks.com/matlabcentral/fileexchange/46261-the-shortest-distance-from-a-point-to-ellipsoid
• BURSA, M; SIMA, Z. Tri-axiality of the Earth, the Moon and Mars. Stud. Geoph. et Geod, 24, (3), 211-217, (1980).
• EBERLY, D. Least Squares Fitting of Data, Geometric Tools, LLC, http://www.geometrictools.com . (2008).
  » http://www.geometrictools.com
• FELTENS, J. Vector method to compute the Cartesian (X, Y , Z) to geodetic (j, λ, h) transformation on a triaxial ellipsoid. J Geod, 83:129-137. (2009).
• LIGAS, M. Cartesian to geodetic coordinates conversion on a triaxial ellipsoid. J. Geod., 86, 249-256, (2012).
• LI ,Q.., GRIFFITHS, J.G.. Least Squares Ellipsoid Specific Fitting. Proceedings of the Geometric Modelling and Processing 2004 (GMP'04). (2004).
• PANOU, G.; DELIKARAOGLOU, D.; AND KORAKITIS, R. Solving the geodesics on the ellipsoid as a boundary value problem. J.Geod. Sci., 3, 40-47, (2013).
• PETROV, Y. (2009): Ellipsoid Fit, 2009. Access on: 07/07/2014. http://www.mathworks.com/matlabcentral/ fileexchange/24693-ellipsoid-fit.http://en.wikipedia.org/wiki/Ellipsoid. Access on: 05/06/ 2014
  » http://www.mathworks.com/matlabcentral/

Publication Dates

• Publication in this collection
  Oct-Dec 2014

History

• Received
  July 2014
• Accepted
  Aug 2014