Unmanned aircraft navigation for shipboard landing using infrared vision

This paper addresses the problem of determining the relative position and orientation of an unmanned air vehicle with respect to a ship using three visible points of known separation. The, images of the points are obtained from an onboard infrared camera. The paper develops a numerical solution to this problem. Both simulation and flight test results are presented.


I. INTRODUCTION
The Brown Water Navy doctrine recently issued by DoD calls for the Naval ships to operate in a close vicinity of the enemy shore. Furthermore, the Navy foresees increasing reliance on the use of unmanned air vehicles (UAVs) for reconnaissance and other missions. This combination of greater utilization of UAVs while operating in proximity of the enemy has highlighted the necessity of using stealth in the recovery of UAVs by the Naval ships. In other words, the ship will not jeopardize its security by communicating its position to the UAV. Therefore, this consideration rules out employing such position sensors as GPS and places emphasis on other passive sensors. Clearly, the only passive sensors capable of providing relative position information are vision based. Moreover, since UAVs are expected to operate around the clock and in all weather conditions infrared (IR) cameras are the passive sensors of choice.
The UAV shipboard autoland includes finding the ship, constructing a landing trajectory based on the relative position, velocity and orientation information (the navigation solution) obtained from passive sensors, and then tracking this trajectory using onboard control system. We are not considering the landing itself because at short ranges low-power communication between UAV and the ship is allowed.
Determining the navigation solution with respect to the ship using passive sensors can be divided into two distinct phases: 1) at large distances the ship is seen as a single hot spot by the onboard IR camera; 2) at closer distances additional features can be determined (see Fig. 1). Phase 1 has been addressed in our previous work, see [1,2], where new filtering algorithms have been developed that integrate IR and inertial navigation system (INS) sensor systems to obtain relative position and velocity of a UAV with respect to a ship. These algorithms are designed to handle out-of-frame events and occlusions that are common to vision sensors. In this work we obtain the navigation solution for the Phase 2. Specifically, the problem of determining range and orientation of an aircraft to a ship, which has a minimum of three identifiable points is addressed. The filtering algorithms developed in [1,2] that use a single point have been used to initialize the algorithm developed here. Once again this algorithm is supposed to help bring the UAV as close to the ship as possible (while all three reference points are visible in the IR camera). After this any final landing procedure may be implemented.
The visual range at which a ship on the surface of the sea can be located depends on the contrast of the ship to its surroundings according to Koshmieder's relation [3,4]. The visual range centered around 0.55 micrometers depends on weather conditions and ambient light, which can produce glare obscuring the ship. Although one can often detect ships at great distances in daylight, the visible spectrum is severely limited under poor weather conditions and especially at night for military vessels without running lights [5]. Since all powered ships will radiate strongly in the 8-12 micrometer IR atmospheric window, it is preferable to locate ships by their hot smokestack and engine. Use of the IR greatly simplifies the problem of locating a ship and reduces susceptibility to glare.  The visible spectrum may be used to supplement the IR for long-range detection in clear daylight conditions.
The hot IR smokestack reliably gives the location of the ship with minimal image processing as compared with visible light. Examples of the information available at different ranges are presented in Fig. 1 in a contour plot and surface plot overlaying the image. At the limit of the detection range for a given minimum resolvable temperature difference, the only detectable point may be the smokestack ( Fig. 1(a)), but at closer range there will be many relatively hot points that can be recovered using image processing ( Fig. 1(b),(c)).
However, only a few of the hottest points projected onto the focal plane of the IR camera can be identified reliably at great distances. This observation naturally leads to the following critical question: what is the minimum number of known reference points (RPs) that is necessary to determine the range and orientation of the IR camera with respect to the ship. It turns out the answer to this question is three. But using only three points always results in more than one solution as has been shown by a number of researchers in the areas ranging from projective geometry to photogrammetry. Indeed, a survey of the scientific literature reveals that the number of possible solutions may range from four to fifteen. (A detailed discussion of these results is given in Section III.) This problem of nonuniqueness is usually resolved at very close ranges by using more than three points that must lie in the same plane (see for example [6][7][8] and references therein). However, at greater ranges the computational cost of obtaining each additional known point that lies in plane defined by the initial three becomes prohibitive for real-time applications. Therefore, we assume here that three reliable points may be computed from the location of the smokestack and the extents (width and height) of the ship. An illustration of how these images of three RPs can be obtained from an IR image of the ship is given next. Fig. 2 presents examples of IR images of the naval ship passing through San Diego harbor at two different distances of more than two miles. The images are shown in contour and surface plot above a false color image to illustrate the information that is available. Using previously developed algorithms [5] the ship may be extracted from the background as shown in Fig. 3. The false-color image of the ship has been automatically located and indicated by a box in the first insert. The next two inserts are the single level binary and gray-scale images of the ship that have been extracted from the background.
Finally, Fig. 4 shows images of three RPs that can be used by the algorithm developed here. These images may be obtained as follows: the smokestack by using thresholding of the image directly, while the other two points by intersecting the images of the edges of the ship's deck.
Having shown how the ship may be located and the three points of information (images of RPs) established, we focus our attention on determining  the range and orientation of an IR camera with respect to the ship using images of three RPs (see Fig. 5). To address the problem of non-uniqueness of the solution we introduce a concept of an admissible solution and using extensive numerical analysis show that under reasonable assumptions on the relative geometry of the camera and RPs there can be at most two such solutions. (The admissible solution implies that the camera is in front of the ship.) Based on this analysis we develop an efficient numerical algorithm that identifies the two admissible solutions and selects the correct one. The utility of the algorithm is illustrated in simulation and using flight test data collected by an IR camera mounted on a small UAV. This paper is organized as follows. Section II contains a mathematical formulation of the problem. Section III describes previous work in this area, which dates back to 1841. Section IV contains numerical analysis of the problem and discusses proposed numerical solution. Results of computer simulation are shown in Section V. Sections VI and VII discuss the flight test setup and present the results obtained using flight test data. The paper ends with conclusions.

II. PROBLEM FORMULATION
Consider Fig. 6. Letp i = fx i , y i , z i g, i = 1,:::,3 denote the vectors connecting the origin of the camera frame O with the three known points P i , i = 1,:::,3. Let d i , i = 1,:::, 3 denote distances between these points: and s i = kp i k, i = 1,:::, 3 denote the norms of the vectorsp i . We utilize the pinhole camera model [9]. Using this model the projection of each RP onto the image plane of the camera with the focal length f has the following form: i = 1,:::,3:  (1) and (2) we obtain nine equations in nine unknowns fx i , y i , z i g, i = 1,:::,3. Using (2) we get By substituting these expressions into (1) and after simple algebra we can reduce (2) and (1) to a set of three nonlinear equations in three unknowns: To simplify notation we rewrite (4) as follows Note, that the coefficients A, B, C,d i , i = 1,:::,3 are strictly positive by construction. Using (5) one can obtain another system of equations better suited for further analysis. First, observe that Now by rewriting system (5) in terms of s i , i = 1,:::,3 we get where cos ® 1 = (p 1 ,p 2 ) kp 1 kkp 2 k , cos® 2 = (p 1 ,p 3 ) kp 1 kkp 3 k , cos ® 3 = (p 2 ,p 3 ) kp 2 kkp 3 k (see Fig. 6).
Geometrically system (7) can be described as an intersection of three orthogonal elliptic cylinders with the semiaxes rotated around corresponding symmetry axes by the angle of 45 ± . This follows directly from the canonical form of equation (7). The magnitudes of the semiaxes for each cylinder are equal to It is clear that the intersection of any two cylinders is always nonempty and the number of solutions in this case is infinite. However, by adding a third cylinder one can get only a finite number of intersection points. In practice for the system (7) this number cannot be zero or two (as will be shown in Section IV). The only possible set of solutions contains four, six, or eight points. For instance, Fig. 7(a) demonstrates an example with four real solutions to system (7) (two pairs of symmetric points). Increasing the size of the cylinder along the s 2 axis, results in three pairs of solutions ( Fig. 7(b)). Further increase leads to four pairs ( Fig. 7(c)) and again to three pairs of solutions ( Fig. 7(d)). Eventually only two pairs of solutions remain ( Fig. 7(e)).
In the work reported here, we make the following assumption.
A1. The camera is always in front of the plane defined by three RPs P i , i = 1,:::,3.
In the sequel the set of all vectorsp i = fx i , y i , z i g, i = 1,:::, 3 that satisfy Assumption A1 is called admissible. Assumption A1 implies that the x-component of each vector in the admissible set is positive (i.e., s i > 0, i = 1,:::,3).
Summarizing, the problem to be addressed here is to find all admissible solutions to (7) using Assumption A1 and two additional assumptions discussed in Section IV. Furthermore, since the set of admissible solutions contains more than one element we would like to develop a test to select the correct solution.

III. PREVIOUS WORK
It turns out that the three-point perspective pose estimation problem (P3P) (as the problem addressed here is called in computer vision) was first formulated by the German mathematician Grunert in 1841 ( [10]). Since then it has been addressed by many scientists throughout the world. As a result it has been well established that the problem does not have an analytical solution and most attempts were directed at getting a numerical one. In the remainder of this section we present a brief overview of the existing results partly based on reference by Haralick, Lee, Ottenberg, and Nölle ( [11], 1991).
According to Müller ([12], 1925) Grunert obtained (7) by simple use of the law of cosines implemented for the corresponding tetrahedron. Most of the work addressing this problem used formulation (7) as well. Grunert himself introduced two new variables to show that with their help system (7) can be reduced to a fourth-order polynomial with respect to v ¤ . Since the coefficients of this polynomial were expressed as complicated functions of the problem data, this polynomial could be neither solved/simplified analytically nor analyzed. This is due to the fact that, in general, the roots of the fourth order polynomial cannot be obtained analytically (see for example [13] and references therein). The same approach was used by Merritt ([14,15], 1949) and independently by Fischler and Bolles ( [6], 1981). With the same substitution and by manipulating different pairs of equations and different multipliers, they reduced the problem to a fourth-order polynomial in terms of u ¤ , rather than v ¤ .
Several attempts to decrease the order of the final polynomial to be solved are known as well ( [11]). For example, Finsterwalder ([16], 1903) instead of finding all roots of a fourth-order polynomial reduced the problem to finding a root of a cubic polynomial and the roots of two quadratic polynomials. Grafarend, Lohse, and Schaffrin ( [17], 1989) applied different transformations to (7) in order to reduce the problem to finding the same roots of a cubic polynomial and the roots of two quadratic polynomials. In [18] Lohse further extended these results to show that admissible solutions can be picked up from as many as 15 solutions provided by his transformations.
Linnainmaa, Hanvood, and Davis ( [19], 1988), using another and more effective substitution, that is reduced system (7) to fourth-order polynomial in s 2 1 . Quan and Lan ( [8], 1999) mentioned that they could use classical Sylvester resultant to rewrite (7) as an eighth-order polynomial in s 1 (fourth-order polynomial in s 2 1 ). Today availability of powerful computers has made it easy to find all possible solutions to either three-or four or eight-order polynomial. Moreover, Haralick, Lee, Ottenberg, and Nölle [11] did compare numerical accuracy of possible solutions obtained using all approaches mentioned in this section (and showed by the way that only approaches by Fischler and Bolles, and Linnainmaa, Hanvood, and Davis involve no singularity in the computation). But the questions of what is the number of admissible solutions for the specific problem geometry and how to select the correct one still have not been completely answered.
To show that at most four admissible solutions can be found, Fischler and Bolles considered a specific case of equilateral triangle P 1 P 2 P 3 (see Fig. 6) where d i = 2 p 3, cos ® i = 5=8, i = 1,:::, 3, i.e., when system (7) becomes singular. For this case they obtained numerically four admissible solutions shown graphically in Fig. 8(a). Note that for case of equilateral triangle, system (7) is reduced to the following By subtracting the second equation from the first, we obtain (s 2 ¡ s 3 )(s 2 + s 3 ¡ 2s 1 cos ®) = 0: Equating s 2 and s 3 in the third equation of system (11) we get Finally, either the first or the second equation provides two solutions for s 1 : Note that due to symmetry, four admissible solutions can be obtained for this particular case (the second multiplier in (12) gives the same nonsymmetric roots). Moreover, four admissible solutions exist for a more general case, when only two of the three equations in (7) are singular, i.e., when triangle P 1 P 2 P 3 is isosceles and camera resides in the symmetry plane (these solutions can be obtained in the same manner as was done in (11)-(14)). Now, the natural question to ask is could it still be the general case?
To answer this question Fischler and Bolles propose the following procedure to obtain four admissible solutions (see Fig. 8(b)). Moving along the line OP 1 (see Fig. 8(b)) and fixing the pair of points fP 2 , P 0 2 g and fP 3 , P 0 3 g corresponding to edges d 1 and d 2 , respectively, one obtains four candidate solutions ([P 2 ; P 3 ], [P 0 2 ; P 0 3 ], [P 0 2 ; P 3 ], and [P 2 ; P 0 3 ]). In the general case, these candidate solutions have different lengths that do not equal d 3 . Therefore, Fischler and Bolles suppose that each of them can be equal to d 3 for a certain s 1 , resulting in four admissible solutions. That is true. However, this is true only for singular cases discussed above and for a certain configuration of problem data in the general case (as shown in the next section).
The P3P is also well known in the field of photogrammetry. The term photogrammetry came into general use in the U.S. about 1934, although the term already had been widely used in Europe since 1893 after German A. Meydenbauer [20,21]. The main objective of photogrammetry is to obtain reliable landscape measurements by means of aerial photographs. Since tilted photographs introduce errors in the map position, it is important to take into account tilt and swing in aerial photographs at the time of exposure. This task is one of the fundamental problems in the photogrammetry and is called "space resection involving the determination of the spatial position of a camera exposure station." Despite numerous attempts to solve system (7) analytically, the only solution that is in use in photogrammetry today is a numerical solution developed by E. Church in the mid 1930s [22,23,9].
Church's approach considers two pyramidsground pyramid with three RPs and exposure center of a camera, and image-plane pyramid with three image points and the same center (see Fig. 6). This procedure finds a solution that makes the two pyramids coincide and can, in fact, be interpreted as a well-known method of Newton's iterations. This procedure seems to work well when the initial guess is sufficiently close to the correct admissible solution. However, Church's method does not address the issue of the nonuniqueness of the solution. It improves on an initial guess, which has to be quite accurate (within a few percentage points of the true solution [24,20]). Otherwise Church's method does not guarantee convergence to the correct solution. Moreover, even with a good initial guess this method converges only if the exposure station is "high over the ground and located inside the cylinder or the sphere containing three RPs" [20].

IV. NUMERICAL ANALYSIS
In this section we present results of the numerical analysis of the P3P. The critical issue addressed includes determining geometry of the feasibility regions (i.e., the regions that have two, three, or four admissible solutions). First these regions are computed for the example used by Fischler and Bolles [6], who considered the case where RPs form an equilateral triangle. Then a similar exercise is performed for the case of an arbitrary triangle. This exercise suggests that the shape of the feasibility regions is complex and invariant of the shape of the triangle formed by the RPs. Finally, an algebraic analysis motivated by the geometric solution proposed by Fischler and Bolles is given. This analysis is used to develop an efficient numerical algorithm to solve the P3P. Fig. 9 contains the results of a numerical analysis of the specific example of an equilateral triangle given in Fischler and Bolles [6]. Here the shaded areas represent the set of points where four admissible solutions (FASS) to the P3P exist. In particular, left-bottom part of Fig. 9 shows the 3D view of this solution set parameterized by the elevation above the plane formed by the triangle. Other plots in the figure show the top view of each level set. Fig. 9 clearly shows that FASS can be obtained not only at the point of complete symmetry (central point on Fig. 9 at z = 0:2, :::,5:0) as shown in [6], but at many other points as well. In fact, one can see that the geometry of FASS is fairly complex. This explains why it has been so difficult to characterize FASS analytically and why Church's iterative procedure does not converge on its boundary (see Section III).
In Fig. 10 FASS is computed for the case where RPs form an arbitrary triangle. It suggests that the shape of FASS is invariant of the shape of the triangle formed by RPs. Together, Figs. 9 and 10 suggest that FASS is a complex inverted "pyramid" that is normal to the plane formed by the three RPs. At its base it forms a circle that contains the triangle generated by the RPs. Now we make the following additional assumptions.  A2. min i=1,:::,3 s i À max i=1,:::,3 d i . A3. The camera resides outside of FASS. Assumption A2 implies that the camera is sufficiently far from the ship. As is shown next, Assumptions A1-A3 guarantee that the P3P addressed in this work has only two admissible solutions.
First observe that inside the region that satisfies A1-A3 the following inequalities hold Now from first two equations in (7) we obtain the following expressions for s 2 and s 3 : From (16) it is clear that the set of all possible admissible solutions for s 1 lies in the following interval Furthermore, since due to (15) sin ® i > 0, i = 1,2, this interval is never empty. We first consider the case By substituting the expressions for s 2 and s 3 given by (16) into the last equation of (7) we obtain four equations in s 1 . Let we obtain admissible solutions for s 1 . Consider Fig. 11, which includes the plots of ¢ ¡+ , ¢ ++ , ¢ +¡ , and ¢ ¡¡ versus s 1 . It can be seen that solving (19) for s 1 will result in two sets of solutions, one of which is admissible.
In Fig. 12 the area that contains two admissible solutions in Fig. 11 is magnified. Clearly, the set of admissible solutions for s 1 may contain one or two elements. One-element case results from the fact that either s 2 or s 3 in (16) have one solution, which leads us to the conclusion that s 1 = s ¤ 1 . On the other hand, in the two-element case both s 2 and s 3 have two solutions. Finally observe that due to Assumption A2 none of the following expressions can be zero when evaluated at s 1 = 0: = 2d 1 d 2 (cos P 3 P 1 P 2 + cos® 3 ) > 0: (20) Thus, the functional dependence of the expressions ¢ ++ (s 1 ), ¢ ¡+ (s 1 ), ¢ +¡ (s 1 ), ¢ ¡¡ (s 1 ) on s 1 will always have the form shown in Figs. 11 and 12, i.e., when Assumptions A1-A3 hold only two admissible solutions can only be obtained.
Next we consider the case Using previous arguments we conclude that This leads to Substituting these expressions into the last equation of (7) we obtain which results in a unique solution for all s i , i = 1, :::, 3. Is this geometrically possible? Notice that the expression (21) shows that the angles between edges s 2 and d 1 , and s 3 and d 2 are 90 ± (see Fig. 6). But since the solution is unique, by expressing s 1 in terms of s 2 and s 3 , due to symmetry we obtain that s 1 = cos® i¡1 s i , d i¡1 = sin® i¡1 s i , i = 2, 3. This leads us to the conclusion that the angles between the edges s 1 and d 1 and s 1 and d 2 are 90 ± as well, i.e., the triangle contains two right angles. Therefore, the case with one admissible solution is not realizable (it means that the camera is located at infinity with respect to the plane generated by the three points). This statement can also be proved algebraically. Replace ® 3 in the right hand side of the denominator in (22) with ® 1 + ® 2 . Because of (15), ® 1 + ® 2 > ® 3 . Using this substitution the denominator can be reduced to sin j® 1 ¡ ® 2 j. Thus, the following set of equations should hold However, for small values of angles ® 1 and ® 2 this leads to the inequality jd 1 ¡ d 2 j¸d 3 , which implies that in this case cos P 3 P 1 P 2¸1 (see Fig. 6). Therefore, the region outside of FASS that satisfies Assumptions A1-A3 contains two admissible solutions. Now for completeness we analyze what happens to the solution of the P3P as the camera traverses FASS. Consider Figs. 13(a)-(b). The plots shown here were obtained for the points A through F defined in Fig. 9 in the same way as the plots shown in Fig. 12 (recall the intersection of each graph with the x-axis determines an admissible solution).
Notice the points A and F lie outside of FASS and result in only two admissible solutions; point B lies on the boundary of FASS and produces three admissible solutions (it corresponds to the cases shown in Figs.  7(b) and 7(d)). The rest of the points lie inside of FASS and result in four admissible solutions.
Since we have shown that under Assumptions A1-A3 the nonsingular system (7) always has two admissible solutions, this result can be used to develop a test that selects the correct solution. The idea is to numerically obtain both admissible solutions to (7). Then using (6) and (3) construct two sets of vectors p i 1 orp i 2 , i = 1,:::, 3. Then, compute the normals to the plane generated by each solution and utilize them to identify the correct one. Now, is it possible that two admissible solutions have the same normal? It is fairly easy to see that both solutions have the same normal if they are colinear, i.e.,p i 1 = ¹p i 2 , i = 1,:::,3, since in this case solutions must lie on parallel planes. Now by applying condition (1) we deduce that ¹´1. Thus two admissible solutions always have different (noncolinear) normals. Therefore, the correct solution can be determined by analyzing normals generated by each solution. Using normals to resolve ambiguity is a standard device employed in the structure from motion literature (see for example [25]). This test will fail to identify the correct solution for the case where camera resides inside or on the boundary of FASS. This can be clearly seen for points B through E in Fig. 13, where two or more of the solutions are very close. This implies that resulting normals will be almost colinear. This observation underscores the importance of Assumption A3 for the algorithm presented next.
Based on the results presented above we propose the following algorithm for solving the P3P. Suppose a good initial guess of the normalñ (0) to the plane generated by the three points is available. Then, for step k: 1) solve numerically (10) for x (k) 1 in the interval (17), using x (k¡1) 1 as an initial guess, 2) substitute each solution x (k) 1 obtained in 1) into (3) to getp Using the solution provided by the P3P algorithm the relative orientation of the aircraft with respect to the plane formed by the three RPs can be computed as follows [26]. Let f3pg denote an orthogonal coordinate system attached to the plane generated by the three RPs, let fcg denote the coordinate system attached to the camera and let c 3p R be the coordinate transformation from f3pg to fcg. Form three orthogonal vectorsr 1 ,r 2 ,r 3 using the correct solutionp 1 ,p 2 ,p 3 as follows: Then c 3p R = [r 1r2r3 ]. The transformation matrix c 3p R can also be expressed using Euler angles c 3p R = cos Ã 3p cos µ 3p sin Ã 3p cos µ 3p ¡ sin µ 3p cos Ã 3p sin µ 3p sin Á 3p ¡ sin Ã 3p cos µ 3p sin Ã 3p sin µ 3p sin Á 3p + cosÃ 3p cos Á 3p cos µ 3p sin Á 3p cos Ã 3p sin µ 3p cos Á 3p + sin Ã 3p sin µ 3p sin Ã 3p sin µ 3p cos Á 3p ¡ cos Ã 3p sin Á 3p cos µ 3p cos Á 3p where Ã 3p , µ 3p , Á 3p are yaw, pitch, and bank angles, respectively, with respect to the plane formed by the three RPs. Therefore, one can easily find the Euler

V. APPLICATION TO SHIPBOARD NAVIGATION
Next we present a simulation example where the P3P algorithm is applied to the determination of the range of the aircraft with respect to the ship.
The simulation scenario is shown in Figs. 14-15. The ship is moving North at a constant speed of 10 m/s. Its motion is characterized by pitch and heave oscillations with a period of 12 s. The aircraft is performing a left turn with descent from the initial point (¡1450, ¡200, 470) m with respect to the ship's initial position at an airspeed of 53 m/s. The camera's focal length is f = 0:1 m and declination angle with respect to aircraft longitudinal axis is ¡6 ± . The errors in the projection of each RP onto the image plane of the camera are modeled as independent Gaussian random process with zero mean and standard deviation of one pixel.    Fig. 15 gives the corresponding 3D representation. Fig. 16 contains elliptical cylinders with the coefficients computed from the data taken at the 10 th second of the simulation. Fig. 17 includes the time histories of the z-components of the normals generated by each solution. The z-component of the normal corresponding to the correct solution is close to ¡1 (in camera coordinate frame). This figure also shows how the algorithm switches between the two solutions. Fig. 18 shows the differences between true and estimated values of the components of the vectors p i = fx i , y i , z i g, i = 1,:::, 3 versus relative range to the ship. Clearly the errors are decreasing as the aircraft approaches the ship.

VI. FLIGHT TEST SETUP AND DATA
Naval Postgraduate School has recently completed the development of a rapid flight test prototyping system (RFTPS) for a prototype UAV named Frog Fig. 18. Errors/range history. [27]. The RFTPS consists of a test bed UAV equipped with an avionics suite necessary for autonomous flight, and a ground station responsible for flight control of the UAV and flight data collection, as shown in Fig. 19. A functional block diagram of the RFTPS is also shown in Fig. 19.
The RFTPS provides the following capabilities. Within the RFTPS environment, one can synthesize, analyze, and simulate guidance, navigation, control, and mission management algorithms using a Fig. 19. RFTPS at Naval Postgraduate School. high-level development language; algorithms are seamlessly moved from the high level design and simulation environment to the real time processor; the RFTPS utilizes industry standard I/O including digital to analog, analog to digital, serial, and pulsewidth modulation capabilities; the RFTPS is portable, easily fitting into a van. In general, testing will occur at fields away from the immediate vicinity of the Naval Postgraduate School; the UAV can be flown manually, autonomously, or using a combination of the two. For instance, automatic control of the lateral axis can be tested while the elevator and throttle are controlled manually; all I/O and internal algorithm variables can be monitored, collected, and analyzed within the RFTPS environment.
To test the developed navigation algorithms the Frog UAV was equipped with an Infrared Components Corporation MB IRES IMAGE CLEAR TM Uncooled Microbolometer Module-based IR camera. The camera included a Boeing U3000A uncooled 8-12 ¹m sensor and the Microbolometer Module which produced National Television Standards Committee (NTSC) video signal and output it via a RS-232 interface. The focal length of the camera lens as installed in the Frog UAV was 25 mm with a field of view of 40 ± £ 30 ± . The pixel resolution of the camera video was 320 £ 240. The camera was shock mounted in the nose of the aircraft (see Fig. 19), and the pointing angle was fixed in the x-z plane of the aircraft body axes, declined 5 ± from the longitudinal axis of the aircraft. As a result of the fixed mounting in the aircraft, the aircraft heading and attitude alone determined the camera pointing angle. Because the focal length of the camera was fixed, the camera's field of view was fixed.
The IR camera video signal was recorded in the Frog using a Sony Digital Video Walkman, model GV-D300. This digital video tape recorder (VTR) recorded the live NTSC format video signal from the IR camera in the Digital Video (DV) Standard format on a DV mini video magnetic cassette using a helical scan. In addition to the video image, the VTR also recorded the elapsed recording time of each video frame. Flight tests in support of this project were conducted at an airfield at Camp Roberts, CA. Three charcoal grills were used to simulate the hot spots on the ship (see Fig. 20).
To determine the precision of the P3P algorithm the Frog UAV was equipped with a Trimble AgGPS  132 Differential Global Positioning System (DGPS). The AgGPS 132 system consisted of a 12 C/A-code channel receiver, a combined GPS/DGPS receiver, and a ruggedized antenna. The receiver included ground beacon and satellite DGPS capability. The receiver produced messages that included aircraft latitude, longitude, antenna height (altitude), GPS quality indication, number of satellites, horizontal dilution of precision, speed over ground, and magnetic variation. These messages were transmitted in ASCII format via 38KBaud spread spectrum radio frequency data modems to the ground station. Samples of UAV trajectories recorded by onboard DGPS are shown in Fig. 21. The data obtained by DGPS together with the WGS-84 coordinates of the charcoal grills was used to evaluate the accuracy of the P3P solution during post flight analysis.

VII. FLIGHT-TEST DATA ANALYSIS
The landing sequence was digitized using a frame grabber at a rate of 30 Hz (see Fig. 22) and a simple image-processing algorithm was developed to identify three hot spots on the runway and implemented in real-time on a Pentium-II PC.
The image processing problem, i.e., that of finding the hot spots in the image on the runway, turned out to be nontrivial due to the presence of multiple hot spots in the surrounding area. This is in contrast to finding hot spots on a ship, where they are clearly much hotter than the ocean (compare plots on Fig.  23).   As a result an image processing algorithm was developed to find and track the hot spots observed by the IR camera onboard the UAV Frog [28]. The algorithm consisted of two steps. The first step included finding the hot spots in the initial image and involved a search over the complete image plane (see Fig. 24). The cornerstone of the first step included 1) computing a running average of each row of pixels, 2) subtracting the average from the actual value of each pixel, and 3) and selecting the point that exceed 10 sigma from the average. Once the hot spots were found in the initial image, they were tracked for the remainder of the approach (see Fig. 25). The tracking algorithm involved 1) computing a bounding box around the hot spots, 2) using thresholding and Gaussian weighting to identify the groups of hot spots corresponding to each image of RP within the bounding box, and 3) using inertial data to predict the approximate location and size of the bounding box in the next image. This data plus the geometry of three hot spots with respect to the runway were used further by the P3P algorithm to determine range and orientation of the Frog with respect to the runway.
Results of this application are shown in Fig. 26 together with the trajectory obtained by the onboard DGPS. Fig. 27 shows the errors between the DGPS positions and the trajectories computed using P3P algorithm. Clearly, the algorithm performed as expected with the total error decreasing as a function of range.

VIII. CONCLUSIONS
In this paper we have developed a new efficient solution to the P3P problem and applied it to UAV navigation relative to the shipboard. We have shown that in this particular application the problem has only two admissible solutions, and we have developed a numerical algorithm that determines both. A simple test was proposed to select the correct admissible solution. This numerical algorithm was tested in simulation and using flight test data. The accuracy of the resulting solution was evaluated using onboard DGPS system. The algorithm was determined to perform well. Finally, the algorithm was implemented on a Pentium II computer and was successfully tested in real-time using the flight test data provided by the onboard IR camera at 20 Hz.