ISSN : 2093-8470(Online)

**Journal of Navigation and Port Research Vol.38 No.2 pp.89-96**

DOI : https://doi.org/10.5394/KINPR.2014.38.2.89

# Development of Computer Program for Solving Astronomical Ship Position Based on Circle of Equal Altitude Equation and SVD-Least Square Algorithm

### Abstract

### 1. Introduction

Celestial navigation has been one of the traditional methods to find out the ship's position in the sea. When ships navigate through open sea where it is hard to find landmarks, this method is used as a unique method in determining the ship's position and the compass error. Nowadays, by the rapid development of marine technology science, modern equipment has been invented to calculate ship position accurately, effectively and quickly such as Global Positioning System (GPS), Inertial Navigation System (INS).ie. However, it does not mean that there is no accidents happening during operation time while cut off power, and failure in signal transceiver are possible. Therefore, celestial navigation has continuously been taught for navigational students, deck officers because it is ruled in STCW code 78/95. In 2010, at the IMO conference organized in Manila, a new amendment regarding Seafarers’ Training, Certification and Watchkeeping (STCW) code was carried out, in which, the celestial navigation was recommended to use with electrical nautical almanac and software of celestial navigation calculation. That can be seen as an effort to continue the research and development of celestial navigation into a modern method. The idea using intercept method for celestial navigation was proposed by Marcq de St. Hilaire in 1875. This method determines parameters of Line of Position (LOP) based on observing the bodies’s altitude, calculating the azimuth and the bodies’s altitude at assumed the ship position, then the LOP is ploted on Mercator chart while the intersection of LOPs gives out the ship position. To eliminate working of the plotting, (Dewit, 1974) used least square problem in multi-star observation. In addition, this method requests some constraint conditions to ensure the ship position’s accuracy as following: co-altitude does not exceed 25 minutes, the distance between the true ship position and the assumed ship position is in range 30 NM. However, as LOP is a tangential straight line which is near or on the circle of equal altitude (COP), so the ship position given LOP has low accuracy. The Sphere Triangle Method (STM) can calculate the ship position directly without an additional graphical procedure, however, all formulas adopted for the steps of the computation procedure are necessary in order to seek for an optimal procedure. Moreover, although the STM is a direct computational method, the solving process is still indirect (Chen, 2003).

Recently, some of new methods have been proposed such as SEEM (Chen, 2003), genetic algorithms (Ming-Cheng, 2012), in which the both of them were applied to directly calculate for astronomical ship position. While the SEEM requests a complicated procedure for determining ship position, genetic algorithms is a new approach to celestial navigation, however, there are still some problems in them.

In the vector-matrix method, solving COP system with least square problem and converting into sphere system to give the ship position. However, this method was only applied to situations which there are more than 2 bodies observations. In addition, in solving over-determined equation system, it is possible to meet some problems such as there is no solution to equation system, celestial matrix is rank deficient and numerical stability is lowest when normal decomposition is used.

In this paper, the authors propose some improvements to enhance the accuracy of ship position based on circle of equal altitude equation on Cartesian coordinate system, Singular Value Decomposition–least square algorithm and virtual celestial body. Therefore, astronomical ship position will give more accuracy than previous methods. Finally, by proposal algorithm, a computer program have been developed to help the navigators in calculating directly the ship position when the modern equipment has failure.

### 2. Establishing circle of equal altitude equation on Cartesian coordinate system

In spherical coordinate system, circle of equal altitude equation can be expressed in the Eq.(1)

Where H_{o}, φ, δ, t_{L} are respectively observation altitude, position's latitude, declination and local hour angle of celestial body.

In navigation, the Eq.(1) represents the exact locus of ship position at time of body observation; however, it is not used because of complicated solving procedure. Consequently, the intercept method applied Taylor’s analysis to find out LOP equation which is an equation of tangental line to COP, as shown in Fig. 1, and the LOP can be seen as approximation of COP. Therefore, the system utilizing LOP for ship position are simpler than that of COP equation.

**Fig. 1 **The COP and LOP for astronomical ship position

By using LOPs, the ship position will be determined with lower accuracy than COPs. In this paper, the authors propose another way to establish circle of equal altitude equation on cartesian coordinate system instead of one on spherical coordinate system. As proposed COP equation, the system for ship position is easy to make the solution. The detail of COP can be illustrated in Fig. 2.

**Fig. 2 **The illustration COP in Cartesian coordinate system

In Fig. 2, in order to simplify this problem, the radius of the celestial sphere is chosen as (R = 1). The XOY, XOZ, C_{i}, P_{s}, A are respectively the equator celestial plane, the greenwich celestial plane, the position of celestial body, the ship position and the intersection of COP's plane and .

The vectors of celestial body and ship position can be written as [x_{i}, y_{i}, z_{i}]^{T} ; [X, Y, Z]^{T}. They are expanded into

Where x_{i}, y_{i}, z_{i}, X, Y, Z, λ, t_{G} are cartesian coordinate of celestial body and ship position, longitude, greenwich hour angle of celestial body. And variables are also ruled as follows^{}

Circle of equal altitude is the intersection of its plane and celestial sphere. The plane’s equation of COP which has normal vector of plane is shown in Eq.(2).

Where (a_{o},b_{o},c_{o}) is cartesian coordinate of A (intersection of and COP's plane), therefore

Where k is scale coefficient of vectors was determined according to Eq.(5)

As the celestial bodies lie on the surface of celestial sphere, their co-ordinates can be expressed as Eq.(6)

Applying the Eq.(4), Eq.(5) and Eq.(6) into Eq.(2), the final equation of COP can be explained as Eq.(8)

The Eq.(8) can be used to calculate ship position directly in Cartesian coordinate system. If the spherical coordinate system are applied into Eq.(8), the Eq.(1) can be obtained.

### 3. The SVD-Least square algorithm for solving astronomical ship position

#### 3.1 Virtual celestial Body

According to vector-matrix method, it is very difficult to implement in the cases of 2 celestial bodies. In this paper, the assumed ship position have been used as a virtual body whose distance to the true ship position is co-altitude, the relation between the virtual celestial body and the ship position can be presented in Fig. 3.

**Fig. 3 **The relation of virtual body and the ship position

The altitude of virtual body is determined as Eq.(10).

Where h_{ν}, φ_{ν}, _{ν}, t, a, b are respectively the altitude of virtual body, the latitude of virtual body, ship speed, observation time period, big and small radius of earth.

Therefore, the additional equation for equation system when the number of celestial observations are less than three has form as Eq.(11).

Where x_{ν}, y_{ν}, z_{ν} are the cartesian coordinate of virtual celestial body.

#### 3.2 SVD–least square algorithm for ship position

Astronomical ship position is the intersection of circles of equal altitude together. Because, there are always errors in observation altitudes, so that the equation system for the ship position can be expressed as Eq.(12).

The Eq.(12) can be abbreviated into Eq.(13) as follows

Where is the celestial body matrix, , , are respectively the vectors of the ship position, the errors and the observation altitudes.

The Eq.(13) is solved by the least square method as Eq.(14)

Apply the normal decomposition in Eq.(15), the solution is shown as Eq.(16)

For solving the overdetermined system Eq.(12), the normal or QR decomposition can be used in A matrix. In this research, the authors propose to apply singular value decomposition-Least square algorithm in order to improve numerical stability, it means that the solution of proposed algorithm will give more accuracy than normal method.

Apply the singular value decomposition to A matrix as presented in Eq.(17). Then, the matrices of Eq.(15) will be expanded appropriately to get the solution.

Where U_{(n×n)} is the columns orthonormal matrix satisfying the condition U^{T}U = I_{(n×n)}, its columns are the eigenvectors of AA^{T}. D_{(n×3)} is the diagonal matrix of singular values with the condition D = D^{T}, (σ_{1 } ≥ σ_{2 } ≥ σ_{3}) are the square roots of the eigenvalues of A^{T}A and V_{(3×3}_{)} is the orthonormal matrix satisfying the condition V^{T} = V^{-1}, its columns are eigenvectors of A^{T}A.

In order to find the U matrix, the eigenvalues β_{n} and the eigenvectors of AA^{T} are defined by Eq.(18), Eq.(19)

From the values of (ν_{n1}, ν_{n2}, …, ν_{nn}), the normalizing process is implemented as follows Eq.(20)

then, the U matrix is determined by Eq.(21)

The computing process of V^{T}, D is similar to U, however the matrix is used to defined the eigenvalues and the eigenvectors as A^{T}A. The A^{T}A value is also decomposed as Eq.(22), then it is replaced into Eq.(15).

The Eq.(24) shows the result received from SVD-least square algorithm. By applying Eq.(24), the ship position on Cartesian coordinate system will be more accurate. After that, from this system, the ship position is converted to spherical coordinate system. Finally, they are taken into the earth geodesic system as shown in Eq.(25).

**Fig. 4 **The ship position on the earth geodesic system

Where e = 0.081819791 (WGS-84)

In order to make clearly the calculation procedure of ship position, the flow diagram for proposed algorithm was presented as Fig. 5.

**Fig. 5 **Flow diagram for proposed algorithm

In here, the observation altitudes are taken by sextant at times UT_{1}, UT_{2}, UT_{3}, the parameters for altitude correction are also needed to calculate the observed altitudes. Simultaneously, cartesian coordinate of stars are determined at UT_{1}, UT_{2}, UT_{3}. Solving the equation system with SVD-least square algorithm and converting to spherical coordinate system, it will get the ship position exactly.

### 4. Development of calculation program

After building flow diagram for calculation process, we developed a calculation program for solving ship position by SVD-least square algorithm. The intercept method was also used to calculate ship position, the true position of ship was pre-determined on chart. And then, the accuracy of methods will be compared to one another. Visual Basic software was also used to programming. The interface of program is shown as Fig. 6.

**Fig. 6 **Interface for calculation program

As shown in Fig. 6, the input data for corrected altitude are entered in part of general data such as: course, speed, index error, eyes height, pressure, temperature. The observation time, stars name are taken into Data of equal altitude circles, the Cartesian coordinate of celestial bodies are shown in a table and the ship position are expressed in the results.

In order to verify effectiveness of proposed algorithm, some experiments were performed with numerical database of sun and stars. In case of using the sun or two stars, since the equation system for ship position is insufficient, the virtual celestial body was used. In other cases, altitude observations are more than 2, the virtual celestial body can be also adopted to improve the accuracy of the system's solution.

#### 4.1. Calculation for ship position by sun

The ship position was calculated by sun observation in order to verify the developed program, the input data for calculating ship position were given in Table 1.

**Table 1. **The input data for sun observation

The result of ship position given by SVD-Least square algorithm in case of sun observation was represented in below interface Fig. 7.

**Fig. 7 **The program's interface fix by Sun

After that, the result of ship position is compared to one given by intercept method and true position on chart as shown in Table 2.

**Table 2. **The results of ship position with methods

The Table 2 shows, in case of sun observation, the ship position determined by SVD-least square algorithm get more accuracy than the one by Intercept method. It means that the developed program have good performance and SVD-least square algorithm is higher accuracy and reliability than Intercept method.

#### 4.2. Calculation for ship position by 2 stars

The performance of the developed program was also verified in case of 2 stars observation, the input data for calculating ship position by 2 stars were given in Table 3.

**Table 3. **The input data for 2 stars observation

The result of ship position in case of 2 stars calculated by SVD-least square algorithm is shown in below interface Fig. 8.

**Fig. 8 **The program's interface fix by two stars

The results of ship position are also compared to one given by intercept method and true position on chart, as presented in Table 4.

**Table 4. **The results of ship position with 2 stars

The Table 4 shows the result in case of 2 stars observation as Rasalhague and Vega, the ship position determined by SVD-least square algorithm have also more accuracy than the one by Intercept method. It means that the developed program have good performance and SVD-least square algorithm is higher accuracy and reliability than Intercept method.

#### 4.3. Calculation for ship position by 3 stars

Finally, the ship position was calculated by 3 stars observation as Rasalhague, Vega and Deneb in order to verify the developed program, the input data for calculating ship position were given in Table 5.

**Table 5. **The input data for 3 stars observation

In this case, the virtual celestial body was also applied to improve the accuracy of the ship position. The ship position calculated by proposed method is shown in the interface Fig. 9.

**Fig. 9 **The program's interface fix by three stars

Table 6 shows the ship position in case of 3 stars observation taken to compare the accuracy between methods.

**Table 6 **The results of ship position with 3 stars

The Table 6 shows also the good result of ship position in case of 3 stars observation as Rasalhague, Vega and Deneb. SVD-least square algorithm in this case have also more accuracy than the one by Intercept method. Therefore, the developed program have good performance and SVD-least square algorithm is higher accuracy and reliability than Intercept method in all of 3 cases.

At the end of this experimental part, some summary of position error are shown in the below Table 7.

**Table 7 **Some summary of position error

In this research, the sun's altitude is taken at two different time, in order to overcome this problem, both methods, the sun's altitude at the first time is corrected to the second time, it means that there are two values of the sun's altitude at the second time.

Practically, It always exists the error triangle (cocked hat) in all methods, the least square method is usually used to find the most confident solution.

The ship postion obtainted by two methods is compared to the true ship position. The results show that SVD-least square algorithm and circle of equal altitude has more accuracy than intercept method (LOP).

However, the ship position given on proposed method still has error. To reduce this error, we can decrease error on observation altitude with sextant.

### 5. Conclusion

The conclusion of this study can be summarized as follows

- The circle of equal altitude equation was established on Cartesian co-ordinate system to calculate directly the ship position.

- The SVD-least square algorithm was proposed as a most numerical stable method in order to solve over-determined system for the ship position.

- The calculation program based on suggested algorithm were built to help the navigator in determining astronomical ship position.

- Some experiments were carried out to verify effectiveness of suggested algorithm. The results showed ship position’ accuracy to be higher than traditional methods in all of three cases.

In this research, the authors performed some improvements to celestial navigation method, it is possible to take one more accurately and modernly. In future, if the altitude of celestial bodies are determined continuously by electrical equipments, this method has the potential to be integrated to electrical chart and the e-navigation system.

Although there are some improvements, the celestial navigation method can not replace for GPS. However, with its advantages such as low cost, hight independence and the developments presented in this research, it can be used as backup system for GPS when this equipment has failure.

### Acknowledgement

This research was a part of the project titled "Development of the Intelligent Technology for Marine Accident Protection & Salvage" funded by the Ministry of Land, Transport and Maritime Affairs, Korea

### Reference

2.Dewit C (1974), "Optimal Estimation of a Multi-Star Fix". Journal of The Institute of Navigation, Vol. 21, No. 4, pp. 320 - 325 (1974). http://www.ion.org/

3.Emmett J. Ientilucci (2003), "Using the Singular Value Decomposition", Report of Chester F.Carlson Center for Imaging Science Rochester Institute of Technology, May 29, 2003, pp. 1-8. http://astro.rit.edu/

4.Georgy Gimel' farb (2005), "Least square algorithm". COMPSCI 369 Computational Science, pp. 1-51.

5.Henning Umland (1997), "A short guide to celestial navigation", pp. 1-92. http://www.fabiopeixoto.com/

6.Kaplan G. H. (1996), "Determining the Position and Motion of a Vessel from Celestial Observations", Journal of The Institute of Navigation, Vol. 42, No. 4, Winter 1995-1996, pp. 633-650. http://aa.usno.navy.mil/

7.Kaplan G. H.(1996), "The Motion of the Observer in Celestial Navigation", Navigator's Newsletter, Issue 51 (Spring 1996), pp. 10-14. http://aa.usno.navy.mil/

8.Ming-Cheng Tsou (2012), "Genetic algorithm for solving celestial navigation fix problems", Polish Maritime Research 3(75) Vol 19, pp. 53-59, 10.2478/v10012-012-0031-5.

9.Vulfovich B., Fogile V (2010), "New Ideas for Celestial Navigation in the Third Millennium", The Journal of Navigation, Vol. 63, issue 02, April 2010, pp. 373- 378.: http://dx.doi.org/10.1017/S0373463309990348