The Derivation of the Equations and Algorithm Used to 

Track Communications Platform Locations in the 

Deployment Module of the Network Connectivity 

Analysis Module (NCAM) Software 

by G. Welles Still 

ARL-TR-5217 June 2010 

Approved for public release; distribution is unlimited. 
NOTICES 

Disclaimers 

The findings in this report are not to be construed as an official Department of the Army position unless 

so designated by other authorized documents. 

Citation of manufacturer’s or trade names does not constitute an official endorsement or approval of the 

use thereof. 

Destroy this report when it is no longer needed. Do not return it to the originator. 
Army Research Laboratory 

Aberdeen Proving Ground, MD 21005-5068 

ARL-TR-5217 June 2010 

The Derivation of the Equations and Algorithm Used to 

Track Communications Platform Locations in the 

Deployment Module of the Network Connectivity 

Analysis Module (NCAM) Software 

G. Welles Still 

Survivability/Lethality Analysis Directorate, ARL 

Approved for public release; distribution is unlimited. 
ii

REPORT DOCUMENTATION PAGE Form Approved 

OMB No. 0704-0188

Public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, gathering 

and maintaining the data needed, and completing and reviewing the collection information. Send comments regarding this burden estimate or any other aspect of this collection of information, 

including suggestions for reducing the burden, to Department of Defense, Washington Headquarters Services, Directorate for Information Operations and Reports (0704-0188), 1215 Jefferson 

Davis Highway, Suite 1204, Arlington, VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person shall be subject to any penalty for failing to 

comply with a collection of information if it does not display a currently valid OMB control number. 

PLEASE DO NOT RETURN YOUR FORM TO THE ABOVE ADDRESS. 

1. REPORT DATE (DD-MM-YYYY)

June 2010 

2. REPORT TYPE 

Final 

3. DATES COVERED (From - To) 

September 2006–August 2009 

4. TITLE AND SUBTITLE 

The Derivation of the Equations and Algorithm Used to Track Communications 

Platform Locations in the Deployment Module of the Network Connectivity 

Analysis Module (NCAM) Software 

5a. CONTRACT NUMBER 

5b. GRANT NUMBER 

5c. PROGRAM ELEMENT NUMBER 

6. AUTHOR(S) 

G. Welles Still 

5d. PROJECT NUMBER 

04E5E1 

5e. TASK NUMBER 

5f. WORK UNIT NUMBER 

7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES) 

U.S. Army Research Laboratory 

ATTN: RDRL-SLE-M 

Aberdeen Proving Ground, MD 21005-5068 

8. PERFORMING ORGANIZATION 

 REPORT NUMBER 

ARL-TR-5217 

9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES) 10. SPONSOR/MONITOR’S ACRONYM(S) 

11. SPONSOR/MONITOR'S REPORT 

 NUMBER(S) 

12. DISTRIBUTION/AVAILABILITY STATEMENT 

Approved for public release; distribution is unlimited. 

13. SUPPLEMENTARY NOTES 

14. ABSTRACT 

The Network Connectivity Analysis Model (NCAM) software predicts the viability of links between moving communications 

nodes. The deployment model tracks the position of the moving nodes as a function of time by having the user outline the 

platform paths with waypoints and defines the time the platform leaves and arrives at each waypoint. It then interpolates the 

position of the waypoint at strategically chosen snapshot times. For verification purposes, the derivations of the equations the 

deployment module uses to do the interpolations are presented. The velocity profile assumed was trapezoidal in shape. This 

shape modeled the velocity closely enough so that platform motion was realistic, yet the input required from the user was 

limited and the equations describing the motion were manageable. The results were consistent with the assumed velocity 

profile. NCAM is being developed jointly by the Missile Defense Branch at Aberdeen Proving Ground, MD, and the 

Communications Electronic Warfare Branch at Fort Monmouth, NJ, of the U.S. Army Research Laboratory’s 

Survivability/Lethality Analysis Directorate. 

15. SUBJECT TERMS

great circle, interpolation, oblate spheroid, snapshot time, sphere, waypoint, WGS84, velocity profile, spherical coordinates, 

trigonometric functions, elevation, latitude, longitude, coordinate transfer, arrival time, departure time, cruise velocity 

16. SECURITY CLASSIFICATION OF: 

17. LIMITATION 

OF ABSTRACT 

UU 

18. NUMBER 

OF PAGES 

82 

19a. NAME OF RESPONSIBLE PERSON

G. Welles Still 

a. REPORT 

Unclassified 

b. ABSTRACT 

Unclassified 

c. THIS PAGE 

Unclassified 

19b. TELEPHONE NUMBER (Include area code) 

410-278-3377 

Standard Form 298 (Rev. 8/98) 

Prescribed by ANSI Std. Z39.18 
iii

Contents 

List of Figures v

Executive Summary vii

1. Introduction 1

1.1 Background .....................................................................................................................2

1.2 Purpose ............................................................................................................................2

1.3 Scope ...............................................................................................................................2

2. Methodology 2

3. Analysis 3

3.1 Spherical Latitude-Longitude Coordinate System Review .............................................3

3.2 Equation of a Great Circle on a Sphere ...........................................................................4

3.3 Establishing 0 and max From Two Arbitrary Points .....................................................8

3.4 Special Considerations ..................................................................................................12

3.4.1 (2,2) on the Equator ........................................................................................12

3.4.2 (2,2) or (1,1) but not Both at the Intersection of the Equator and the 

International Date Line or the Prime Meridian .............................................................13

3.4.3 Both (2,2) and (1,1) on the Equator .............................................................14

3.4.4 Both (2,2) and (1,1) on the Same or Opposite Meridians ............................14

3.4.5 (2,2) or (1,1) at a Pole ..................................................................................16

3.4.6 (2,2) and (1,1) Are Antipodal ......................................................................16

3.4.7 (2,2) and (1,1) at Opposite Poles .................................................................18

3.4.8 (2,2) and (1,1) Nearly Antipodal on the Equator .........................................18

3.4.9 Summary of Waypoint Placement Problems .....................................................26

3.5 Value of the Earth’s Average Radius ............................................................................26

4. Results 28

4.1 Definition of Moving Platform Variables .....................................................................29

4.1.1 Mathematical Description of Platform Motion .................................................31

4.1.2 Mathematical Relationship Between the Motion-Defining Variables ..............32
iv

4.1.3 Restrictions on Input Data – General ................................................................34

4.1.4 Distance Calculation ..........................................................................................36

4.1.5 Restrictions on Input Data – Maximum Speed, Minimum Time ......................41

4.1.6 Restrictions on Input Data – Minimum Speed, Maximum Time ......................42

4.1.7 Restrictions on Input Data – Peak Speed Less Than vcru. .................................43

4.1.8 Restrictions on Input Data – Insufficient Time for Platform to Reach Zero .....44

4.1.9 Inputting Speed or Time Data ...........................................................................45

4.1.10 Summary of Inputting Waypoint Data ..............................................................46

4.2 Modifying Input ............................................................................................................49

4.2.1 Moving a Waypoint ...........................................................................................49

4.2.2 Adding a Waypoint ...........................................................................................51

4.2.3 Removing a Waypoint .......................................................................................53

4.2.4 Changing Waypoint Data ..................................................................................55

4.2.5 If a Platform Has Only One Waypoint ..............................................................58

4.3 The Interpolation Process ..............................................................................................59

4.3.1 Calculate the Interpolated Distance ...................................................................61

4.3.2 Calculate the Interpolated Elevation .................................................................61

4.3.3 Calculate the Interpolated Angle i ................................................................62

4.3.4 Calculate the Interpolated Latitude and Longitude ...........................................63

4.3.5 Summary of Inputting Waypoint Data and Obtaining Interpolated Points .......72

5. Conclusions 73

6. References 75

List of Symbols, Abbreviations, and Acronyms 79 

Distribution List 80
v

List of Figures 

Figure 1. Modules of the connectivity software NCAM. ...............................................................1

Figure 2. The latitude-longitude coordinate system superimposed on the rectilinear 

coordinate system.......................................................................................................................3

Figure 3. A sphere with a primed rectilinear coordinate system and great circle added. ...............5

Figure 4. The angle between the great circle and the equator and the planes containing them. 

The sphere has been removed for clarity. ..................................................................................6

Figure 5. Viewing the sphere with planes 1 and 2 of figure 4 along the y’ axis. The sphere 

has been removed for clarity. .....................................................................................................6

Figure 6. Figure 4 with an unprimed rectilinear coordinate system added. ....................................7

Figure 7. View of the geometry in figure 6 in the negative z direction. .........................................7

Figure 8. The function tan  repeats itself every 180°. ...................................................................9

Figure 9. A unit circle showing the change of cos  when  changes by 180°. ...........................10

Figure 10. The two solutions for max and 0 are opposite each other on a great circle. ..............11

Figure 11. The points (max,0), the origin, and (– max0 ± 180°) share the same plane with 

the x’ and the z/z’ axes. ............................................................................................................11

Figure 12. The points where the great circle cross the equator differ in longitude from 

(max,0) by 90°. .......................................................................................................................12

Figure 13. 2 and 2 are the intersection of the international date line and the equator. ..............14

Figure 14. Both (2,2) and (1,1) are on the equator. .................................................................15

Figure 15. The two waypoints are on opposite meridians. ...........................................................15

Figure 16. Two antipodal points do not define a unique great circle. ..........................................17

Figure 17. The shortest path on the oblate spheroid connecting two equatorial waypoints 

nearly opposite each other goes closer to the poles than the equator. The plane provides 

the cross section of the oblate spheroid. ..................................................................................19

Figure 18. The y-z cross section of the oblate spheroid of figure 17. ...........................................20

Figure 19. Cross-sectional view of the equator of figure 15 looking in the negative z 

direction when y0 is at its maximum value and the definition of the angles min and max

added. .......................................................................................................................................21

Figure 20. Two nearly antipodal waypoints on the equator bisected by a plane that is 

perpendicular to the straight line connecting the nearly antipodal waypoints. ........................22

Figure 21. A more general view of figure 17: the straight line connecting the waypoints is 

not necessarily parallel to the x or y axis, yet it is perpendicular to the z axis. .......................23

Figure 22. A more general cross section of figure 18. ..................................................................24
vi

Figure 23. Variables that describe moving platforms as they travel between waypoints. ............29 

Figure 24. Motion profile for moving platform between waypoints. ...........................................31

Figure 25. Two waypoints and the Earth’s center point form a wedge. .......................................36

Figure 26. The wedge from figure 25 with the Earth removed for clarity. ...................................37

Figure 27. The motion profile for maximum arrival speed and minimum arrival time. ...............41

Figure 28. The maximum speed if the distance does not allow enough time for Vmax to equal 

Vcru. ..........................................................................................................................................42

Figure 29. The motion profile for zero arrival speed and maximum arrival time. .......................43

Figure 30. The speed of the platform if it never reaches vcru. .......................................................43

Figure 31. The minimum speed if the platform does not have time to reach zero speed at the 

arrival waypoint. ......................................................................................................................44

Figure 32. The quantity vtop and tt are defined graphically. ..........................................................45

Figure 33. Overview of the waypoint data input process. ............................................................47

Figure 34. Overview of the process of moving a waypoint. .........................................................50

Figure 35. Overview of the process of adding a waypoint. ..........................................................52

Figure 36. Overview of the process of removing a waypoint. ......................................................54

Figure 37. Overview of the process of changing the waypoint data. ............................................56

Figure 38. The process of inputting data for a platform with one waypoint. ...............................58

Figure 39. The process of inputting data for changing the location of a single waypoint 

platform. ...................................................................................................................................59

Figure 40. The process of adding a waypoint to a platform that has one waypoint. ....................60

Figure 41. The first step in transforming the interpolated point to a new coordinate system. .....64

Figure 42. The first step in transforming the interpolated point to a new coordinate system 

with the Earth removed viewing in the –z direction. ...............................................................64

Figure 43. The second step in transforming the interpolated point to a new coordinate 

system: the y’ axis acts as the axel so that y’ = y”. .................................................................65

Figure 44. The second step in transforming the interpolated point to a new coordinate system 

with the Earth removed viewing in the +y’ direction. .............................................................66

Figure 45. The dots on the sine and cosine graph at 0.5 show that there are two distinct 

angles for any one value of sine and cosine in the interval –180° to +180°. ...........................68

Figure 46. For an angle  ranging from –180° to +180°, the sine and cosine functions take on 

either a positive or negative sign in each 90° quadrant. ..........................................................69
vii

Executive Summary 

The Network Analysis Connectivity Model (NCAM) is a software package designed to predict 

the viability of communication links between moving nodes.1

 It consists of six modules: the 

deployment module, the propagation module, the antenna module, the noise module, the link 

budget module, and the connectivity confidence interval (CCI) module. The deployment module 

keeps track of the platform latitude, longitude, and elevation as a function of time. The 

propagation module uses the Terrain Integrated Rough Earth Model2 program to account for the 

effect of terrain on signal propagation. The antenna module applies antenna properties, such as 

gain and radiation pattern, to signal propagation. The noise module accounts for the effect of 

radio noise from man-made sources, such as jamming, transformers, vehicles, transmission from 

other radios, or heat from the broadcasting radio, and natural sources, such as the background 

cosmic radiation, the Sun, and the Earth. The link budget module takes the data from the 

previous modules and calculates the wireless link signal-to-noise ratio for each communications 

link. The CCI module calculates the appropriate standard deviation of the path loss associated 

with each link and computes the probability of successful signal reception. NCAM uses the six 

modules to do the calculations necessary to perform the wireless network simulation. NCAM 

requires all six modules to perform its calculations. 

The approach used is to do all the calculations needed to determine the CCI at certain key times 

called snapshot times. The computations at the snapshot times form a series of samples through 

time of link viability, which, when viewed in sequence, shows how well the nodes are linked 

through time, much like frames in a movie show samples of cinema action. To predict the 

viability of communication links between moving nodes through time, it is necessary to track the 

location of all nodes. 

The deployment module functions by requiring the user to input an outline of the platform’s path 

using a set of sample points called waypoints. Each waypoint is specified with the latitude, 

longitude, elevation (read from a Digital Terrain Elevation Data database for a ground 

platform3), and arrival and departure time of each waypoint. The deployment module then 

interpolates the position of the platform at the snapshot times.

 1Still, G. W.; Nealon, J. F. The Case for Using the Spherical Model to Calculate the Interpolated Points in the Connectivity 

Software Deployment Module; ARL-TR-4373; U.S. Army Research Laboratory: Aberdeen Proving Ground, MD, February 

2008. 

2

Eppink, D.; Keubler, W.; TIREM/SEM Handbook; ECAC-HDBK-93-076; Department of Defense Electromagnetic 

Compatibility Analysis Center: Annapolis, MD, March 1994; pp 1–134. 

3

Pablo’s Mission Planning. http://www.mission-planning.com/DTED_Part2.htm (accessed 28 June 2007). 
viii

For verification purposes, this report documents the derivation of the equations that the 

deployment module uses to interpolate the positions of the moving communications platforms 

(or nodes) and the development of the algorithm for adding, removing, and modifying the 

waypoints. In developing the equations, it was assumed that there were two constant velocities: 

either zero or some maximum cruising velocity vcru. vcru is taken from a database, and its value 

is selected based on terrain roughness and the platform type. When the velocity changes, it 

increases or decreases at a constant acceleration rate of ±a0. This results in a time vs. velocity 

graph that is trapezoidal in shape. Furthermore, the velocity of the platform as it arrives and 

leaves the waypoint must be the same. A platform is allowed to change acceleration from 

+ (or –) a0 to – (or +) a0 as it passes through a waypoint or remain at zero if it is traveling at vcru. 

The platform also has the option of lingering at a waypoint. 

The user is given the choice of specifying platform speed or arrival time at each waypoint. The 

motion equations derived for interpolation also provide a lower and upper bound of the node’s 

speed or arrival time at each platform. So, the user is given the maximum choice for specifying 

node motion but is prevented from entering times and speeds that are impossible for the platform 

to attain. 

The algorithm developed provides for the removal or addition of waypoints, changing the arrival 

time/velocity at waypoints, or changing the location of the waypoints. When the waypoint 

properties are modified, where possible, the waypoint velocity remains the same or is modified 

to be the minimum or maximum allowed for the platform, and the arrival/departure times are 

recalculated for the subsequent waypoints. 

The velocity profile was chosen so as to make the platform motion realistic. A simple profile 

would have been unrealistic. Making the profile more complicated would make the motion more 

realistic, but the derivation of the equations needed for interpolation would have been more 

involved than those presented here. A more complex velocity profile would also have required 

more user input. The velocity profile (and the subsequent derivation of the equations of motion) 

here provided the best compromise between user control, equation and algorithm simplicity, and 

realistic platform motion. 

NCAM is being developed jointly by the Missile Defense Branch at Aberdeen Proving Ground, 

MD, and the Communications Electronic Warfare Branch at Fort Monmouth, NJ, of the U.S. 

Army Research Laboratory’s Survivability/Lethality Analysis Directorate. For verification 

purposes, the details of the algorithm and formulae development have been documented and 

shown to be consistent with the assumptions of the model’s description of platform motion.
1

1. Introduction 

A new software package called the Network Connectivity Analysis Model (NCAM) (1, pp 3–4) 

is currently undergoing development. It is being developed jointly by the Missile Defense 

Branch at Aberdeen Proving Ground, MD, and the Communications Electronic Warfare Branch 

at Fort Monmouth, NJ, of the U.S. Army Research Laboratory’s Survivability/Lethality Analysis 

Directorate. NCAM predicts the viability of communications links between moving nodes. 

Six modules constitute NCAM, as illustrated by the screen capture of the graphical user interface 

shown in figure 1. The deployment module calculates and tracks the movements of the wireless 

network nodes. The propagation module calls on the program TIREM (the Terrain Integrated 

Rough Earth Model) (2) to calculate signal attenuation due to the effect of atmosphere and 

terrain. The antenna module contains data concerning the antennae used, including gain, loss, 

and electromagnetic radiation pattern. The noise module computes receiver noise attributable to 

internal, external, man-made, natural, hostile, and nonhostile noise sources. The link budget 

module takes the data from the previous modules and calculates the wireless link signal-to-noise 

(S/N) ratio for each communications link. The connectivity confidence interval (CCI) module 

calculates the appropriate standard deviation of the path loss associated with each link and 

computes the probability of successful signal reception. NCAM uses the six modules to do the 

calculations necessary to perform the wireless network simulation. 

Figure 1. Modules of the connectivity software NCAM.
2

1.1 Background 

The deployment module was the first module developed. At certain times in the simulation 

called snapshot times, the S/N ratio and the CCI is calculated for every link between the moving 

nodes. To do those calculations, the distance between nodes must be derived from their 

positions. To obtain those distances, the Earth could have been modeled either as a perfect 

sphere (with elevations in the radial direction, whose values are very small compared to the 

radius of the Earth) or as an oblate spheroid (also with small variations in the radial direction to 

represent elevation). Developers chose the spherical model of the Earth as the basis model for 

distance computation because it allowed more rapid computation, with a small, but acceptable, 

sacrifice in accuracy (1, pp 3–4). 

1.2 Purpose 

The user defines the path for each communications node or platform by selecting key locations 

along the path called waypoints. Successive waypoints in time are then connected by great 

circles, which are the shortest lines between two points on the surface of a sphere. The great 

circle approach accounts for the locations and the elevation of each waypoint. 

This report traces the derivation of equations describing platform motion. The equations derived 

will be used by NCAM’s deployment module to obtain platform position as a function of time. 

This allows NCAM to calculate and use the distance information in its calculation of link 

viability. An algorithm allowing for the modification of the number of waypoints for each 

platform and the modification of waypoint location is also included to allow modification of the 

platform paths. This presentation of the motion equations and the algorithm is intended to 

provide evidence for validating the NCAM software. 

1.3 Scope 

The derivations are limited to the deployment module. Because elevations are considered, the 

equations derived can be applied to either ground or airborne nodes. Due to the nature of the 

trigonometric functions, special limitations are also considered with the derivation of the 

equations. 

2. Methodology 

Platform motion is assumed to be along a great circle path. Platform velocities vary between a 

zero or nonzero cruise velocity. Velocities between zero and the cruise velocity are attained 

while accelerating between the cruise and zero velocity at either a positive or negative rate of 

constant acceleration. Solutions to special problems occasioned by the limitations of the 

trigonometric functions are also considered, as is the process needed to add, remove, or modify 

waypoint locations.
3

3. Analysis 

The analysis consists of deriving the equations of motion for the platforms using spherical 

coordinates. After reviewing the spherical coordinate system, the equation for a great circle is 

derived. The equation is derived for two arbitrary points that define the circle, and special 

circumstances such as the two defining points (or waypoints) nearly antipodal on the equator, at 

opposite poles, antipodal, on the same or opposite meridians, one waypoint at the pole, both 

waypoints on the equator, or one point at the intersection of the equator, and the international 

date line, or prime meridian. The analysis ends with the derivation of the value for the Earth’s 

average radius used in the spherical model. 

3.1 Spherical Latitude-Longitude Coordinate System Review 

The latitude-longitude system of coordinates used in modern navigation on Earth is illustrated in 

figure 2.  is the latitude, which varies from –90° to 90°. The designation “north latitude” 

indicates a positive value, while “south latitude” indicates a negative value.  is the longitude, 

varying from –180° to +180°. “East” corresponds to positive latitude, while “west” connotes 

negative latitudes. R is the distance from the Earth’s center and is the sum of the mean sea level 

and the elevation. R is always positive and can have a value between 0 and infinity. 

Figure 2. The latitude-longitude coordinate system superimposed on 

the rectilinear coordinate system.
4

To facilitate angle and distance calculations, we superimpose a three-dimensional (3-D), 

rectilinear coordinate system on the latitude and longitude system. The x axis starts at the 

Earth’s center point and passes through the point where the prime meridian (0° longitude) and 

equator (0° latitude) intersect. Likewise, the y axis starts at the center point but goes through the 

point at the intersection of the equator (0° latitude) and 90° east longitude. The z axis begins at 

the Earth’s center as well but passes through the North Pole (90° north latitude). 

Converting from the latitude-longitude coordinate system to the rectilinear system involves 

employing simple geometric relationships. In the vertical right triangle in figure 2, the side 

opposite the latitude angle  is the z coordinate. Expressed in terms of the latitude-longitude 

system, this becomes 

 z = R sin  (1) 

The side adjacent to the latitude angle  becomes the hypotenuse for the two right triangles in the 

x-y plane shown in figure 2, with a value of R cos . The side opposite the longitude angle  is 

the y coordinate, which, in terms of R, , and  is 

 y = R cos  sin  (2) 

The side adjacent to the latitude angle  is the x coordinate, rendered 

 x = R cos  cos  (3) 

3.2 Equation of a Great Circle on a Sphere 

We now derive the equation of a great circle on a sphere in terms of latitude angle  and the 

longitude angle . First, consider the situation in figure 3: a sphere of arbitrary radius R0 with 

an equator and an arbitrary great circle. Note that the hidden lines are not shown for clarity. A 

great circle is a circle of maximum radius drawn on a sphere so that it too has a radius of R0, 

dividing the sphere into equal hemispheres (3). Note that the equator is also a great circle, and 

that it, the great circle, and the sphere share the same center point. Included is a rectilinear 

primed coordinate system with three axes: x’, y’, and z’. The relationship between the 

nonprimed coordinate system and the primed coordinate system will be established later. The y’ 

axis begins at the center point and goes through one of the two points of intersection between the 

equator and the great circle. The x’ axis begins at the center point and goes through a point on 

the equator such that the x’ axis forms a 90° angle with the y’ axis. The z’ axis begins at the 

center point and goes through the pole, forming a right angle with both the x’ and the y’ axes. 

The sense of the three primed axis is such that x’ cross y’ equals z’. 
5

Figure 3. A sphere with a primed rectilinear coordinate system and great circle added. 

Next, we consider the angle the great circle forms with the equator. Figure 4 shows the great 

circle and the equator with the two intersecting planes: plane 1 contains the equator and the 

sphere’s center point, while plane 2 contains the great circle and the sphere’s center point. (The 

sphere has been removed from figure 4 for clarity.) Planes 1 and 2 form the angle max with each 

other. Plane 2 also forms the angle max with the x’ axis. The x’ axis is also in plane 1, while the 

y’ axis is in the planes’ line of intersection. The z’ axis is normal to plane 1. For the next part of 

the derivation, we consider the appearance of the planes as viewed in the direction of the y’ axis, 

along the line of intersection of planes 1 and 2. 

Figure 5 shows the planes from this viewpoint. We wish to concentrate on the relationship of 

planes 1 and 2 and realize the quantitative relationship between them and the z’ and x’ axes. 

From the point of view of figure 5, plane 2 is a straight line in the x’-z’ axis. So, in the x’-z’ 

system, the line must obey the equation z’ = m x’ + b (4). This is true for every point in plane 2, 

which contains the great circle. Since plane 2 intersects plane 1 at the point x’ = z’ = 0, b 

becomes zero. Because m is the slope of plane 2, this equals the tangent of the angle max. When 

the substitutions established thus far have been made, the equation for a straight line becomes 

 z’ = tan max x’ . (4)
6

Figure 4. The angle between the great circle and the equator and the planes containing them. The sphere has 

been removed for clarity. 

Figure 5. Viewing the sphere with planes 1 and 2 of figure 4 

along the y’ axis. The sphere has been removed for 

clarity. 

Limiting our consideration of the equation for plane 2 to the great circle portion of plane 2, we 

substitute equations 1 and 3 for z’ and x’ and use the value of R0 for the sphere’s and circle’s 

radii in planes 1 and 2. Because we are in the primed coordinate system and the notation must be 

kept consistent, the great circle equation becomes 

 R0 sin ’ = R0 tan max cos ’ cos ’ . (5) 

Equation 5 is the equation for a great circle but only when a great circle intersects the equator in 

the y’ axis. To obtain the more general equation, we superimpose a second set of axes on the 

primed system: the unprimed x, y, and z axes (figure 6). The new z axis is identical in 

magnitude and direction to the z’ axis. Both the z and z’ axes have the same origin point. That 

means that z = z’. Hence, equation 1 shows that the latitude of the primed and unprimed systems
7

Figure 6. Figure 4 with an unprimed rectilinear coordinate system added. 

are the same. In other words,  = ’. Since all the primed and unprimed axes have the same 

origin point, this puts the new x and y axes in the same plane as the old x’ and y’ axes (the 

equator’s plane). The x’ and y’ axes are rotated from the x and y axes by the angle 0. To obtain 

the equation in the unprimed coordinate system, we will view the situation in the negative z (or 

negative z’, since they are the same) direction (see figure 7). 

Figure 7. View of the geometry in figure 6 

in the negative z direction. 
8

Note that plane 2 and the great circle have been removed for clarity. The x (and y) axis makes 

an angle of 0 with the x’ (and y’) axis. An arbitrary angle ’ in the primed coordinate system 

(or  in the unprimed coordinate system) has been added. As is evident from the figure,  in the 

unprimed coordinate system equals ’ + 0. Solving for ’ and inserting the value into equation 

5, taking advantage of the fact that ’ =  and dividing equation 5 by R0 cos ’, we realize the 

equation for a great circle: 

 tan  = tan max cos ( – 0) . (6) 

The point on the great circle that is farthest north of the equator has a longitude of ’ = 0°. This 

follows from the fact that the intersection of the equator and the great circle happens on the y’ 

axis, at ’ = 90° and ’ = –90°. Since the angle between the planes that contains the great circle 

and equator is max, it follows that this is the latitude of the point on the great circle at ’ = 0. 

The fact that this is the most northern point on the great circle and, hence, the point with the 

greatest latitude is the reason it has been designated max in equation 6. As is apparent from 

equation 6, the northern most point on the great circle with latitude max has a longitude of 0 in 

the unprimed coordinate system. 

3.3 Establishing 0 and max From Two Arbitrary Points 

Suppose we are given two arbitrary waypoints with latitude and longitude coordinates (1,1) and 

(2,2). It is possible for them to define a great circle on the sphere. This is because the two 

points and the sphere’s center point define a plane. The same plane intersecting the sphere 

defines a great circle. (Figure 6 shows how a great circle is contained in a plane. This is the 

same great circle shown in figure 3.) 

It is possible to derive the equation for the great circle that connects the two points. By using 

equation 6 with (1,1) and (2,2), we can obtain the values of max and 0. The first step is to 

plug in (1,1) and (2,2) to equation 6, realizing two instances of equation 6 in equations 7 and 

8: 

 tan 1 = tan max cos (1 – 0) (7) 

and 

 tan 2 = tan max cos (2 – 0) . (8) 

We now have two equations (7 and 8) with two unknowns: max and 0. Dividing equation 7 by 

equation 8 results in an expression that has one of the two unknowns, max (and tan max), 

eliminated: 

 tan 1 / tan 2 = cos (1 – 0) / cos (2 – 0) . (9) 

Next, we use the well-known trigonometric identity cos (A – B) = cos A cos B + sin A sin B (5). 

After applying it to the numerator and denominator of equation 9, dividing that expression by cos 

0, then isolating 0 and taking the inverse tangent (arc tan), we find the value of 0:
9

 0 = arc tan [(tan 1 cos 2 – tan 2 cos 1) / (tan 2 sin 1 – tan 1 sin 2)] . (10) 

There are two solutions for 0 such that –180° < 0 < 180° because the function tan  is periodic. 

That is, the value of tan  repeats itself every 180° (6) as illustrated in figure 8. Because the 

possible range of longitude goes from –180° to +180° (i.e., it spans 360°), it is possible to have 

two full periods of the function tan  represented in longitude. Figure 8 shows this. As a result, 

it is possible to have two values of  for any single value of tan . The primary value of 0 as 

expressed in equation 10 will be in either the first or fourth quadrants. This is generally the 

default value given for the arc tan function for most computational software packages (7). It will 

therefore be necessary to add a line of code to make sure that both possible values of 0 are 

found. To find the second value, we subtract 180° if equation 10 gives us a value of 0 that is in 

the first quadrant (i.e., between 0° and 90°) or add 180° if equation 10 produces a value of 0 that 

is in the fourth quadrant (i.e., between 0° and –90°). 

Figure 8. The function tan  repeats itself every 180°. 

0 (along with 2 and 2) will be used to determine the value of max. To find max, we again use 

equation 6. Substituting in the values just derived for 0 with 2 and 2, dividing by cos 

(2 – 0), and taking the arc tan, we find 

 max = arc tan [tan 2 / cos (2 – 0)] . (11) 

Limiting the range of latitude from –90° to +90°, equation 11 gives one value of max for any one 

value of 0. Because our use of equation 10, and adding or subtracting 180° as required, 

produces two values for 0, it follows that there are two sets of values of 0 and max that satisfy 

equations 10 and 11.
10

Cos  is periodic in 360° (8). To see how having two values of 0 affects the values of max in 

equation 11, we consider the unit circle in figure 9. The term cos (2 – 0) suggests that we 

examine how changing the value of 0 by 180° (and hence changing the value of the equation 11 

term 2 – 0 by 180°) changes the value of max. In figure 9, we see that the red line, orange 

lines, and black axis line form two identical right triangles. So, magnitude of the value of cos (

– 180°) for any value of  is the same as the cos . However, the side of the triangle 

representing the cos ( – 180°) is on the opposite side of the origin as cos . Hence, the sign (but 

not the magnitude) of the cosine function changes from negative to positive or positive to 

negative when the argument  changes by 180°. 

Figure 9. A unit circle showing the change of cos 

when  changes by 180°. 

Therefore, changing the value of 0 by 180° in equation 11 changes the sign of the term cos (2 – 0) 

but not the magnitude. Changing the sign, but not the magnitude, of the term cos (2 – 0) 

changes the sign but not the magnitude of max. So, the two sets of values for 0 and max can be 

expressed as (max,0) and (–max0 ± 180°). 

For any point on a great circle with the coordinates (,), the point exactly opposite it is 180° 

away in longitude , equidistant from the equator in latitude  but in the opposite hemisphere so 

it has the opposite sign. The two points opposite each other on a great circle can be expressed as 

(,) and (–, ± 180°). The points (max,0) and (– max0 ± 180°) are two points opposite each 

other on the same great circle. So even though there are two sets of solutions to equations 10 and 

11 for max and 0, they both define the same unique great circle when applied in equation 6. 

This is because of the symmetry of the location the two sets of points on the great circle. 
11

Figure 10 illustrates this symmetry. Note that part of plane 1 has been removed so that the point 

(–max,0 ± 180°) is visible. It is therefore possible for any unique great circle described by 

equation 6 to have two sets of values for 0 and max. To ensure that max is indeed the maximum 

latitude on the great circle, we will choose the value of 0 that makes max positive. Furthermore, 

note that the points (max,0) and (– max0 ± 180°), the z/z’ axis, and the x’ axis, occupy the 

same plane. Figure 11 shows the relationship between the points (max,0) and (– max

0 ± 180°), and the x’and z/z’ axes. 

Figure 10. The two solutions for max and 0 are opposite each other on a great circle. 

Figure 11. The points (max,0), the origin, and (– max0 ± 180°) share the same plane with the x’ and the 

z/z’ axes.
12

3.4 Special Considerations 

Equation 6 is the equation for a great circle in terms of the latitude  and longitude . Equations 

10 and 11 calculate the constants in equation 6 from two arbitrary points on the circle, (1,1) 

and (2,2). There are some values of (1,1) and (2,2) that do not allow equations 10 and 11 to 

produce values for max and 0. We will now review those values. 

3.4.1 (2,2) on the Equator 

This means that 2 = 0. This reduces equation 10 to 

 0 = arc tan ( –1 / tan 2 ) . (12) 

To satisfy equation 12, 2 and 0 must differ by ±90°. Note that the great circle intersects the 

equator at two points. One is at (2 = 0,2), as shown in figure 12. The other is at (2 = 0, 

2 ± 180°). (The second point is obscured by the x’-z/z’ plane in figure 12.) The point (max,0) 

has the maximum latitude on the great circle. The two points where the great circle and equator 

meet have latitude 0°. So, by symmetry, the point with the maximum latitude (max,0) must be 

equidistant from the two points where the great circle and equator meet. That means 

0 = 2 ± 90° and that (max,0) must be in the x’-z/z’ plane, consistent with equation 12. 

Figure 12. The points where the great circle cross the equator differ in longitude from (max,0) by 90°. 

To put it quantitatively and a little more generally, 

 | 0 – 2 | = 90°, or 270° . (13) 
13

Due to the periodic nature of the tangent function (see figure 8), there are two values of 0 that 

satisfy equation 12 in the –180° to +180° range. Note, too, that both values of 0 satisfy equation 

13. Normally, we would use equation 11 and the two values of 0 to find the two corresponding 

values max and take the set of max and 0 that renders a positive max. Unfortunately, cos 

(0 – 2) = 0 because 0 – 2 = ±90°, or ±270° by equation 13. This causes the denominator of 

equation 11 to be zero. Since (2,2) is on the equator, 2 = tan 2 = 0. Thus, equation 11 reduces 

to 0 divided by 0, rendering the value of max as undefined. 

Fortunately, (1,1) is not on the equator. Equation 11 could have been derived to express max in 

terms of 1 by rendering it 

 max = arc tan [tan 1 / cos (1 – 0)] . (14) 

Because (1,1) is not on the equator, 1 is not zero, and consequently, 0 – 1 does not equal 

±90° or ±270°. This fact permits us to use equation 14 instead of equation 11 to calculate the 

two values of max and choose the positive value to make the set (max,0). 

Had it been (1,1) on the equator instead of (2,2), equation 12 would have rendered 0 as 

 0 = arc tan ( – 1 / tan 1 ) . (15) 

Since equation 11 makes no use of 1 or 1, it could have been used to calculate the values of 

max. Hence, in the instance of (1,1) sitting on the equator and (2,2) not sitting on the equator, 

the calculation for 0 and max could have been handled without special consideration. 

3.4.2 (2,2) or (1,1) but not Both at the Intersection of the Equator and the International 

Date Line or the Prime Meridian 

Should (2,2) be at the intersection of the equator and the prime meridian or the equator and the 

international date line, it would not be possible to calculate the value of 0 using equation 10. 

Equation 10 is repeated here as equation 16: 

 0 = arc tan [(tan 1 cos 2 – tan 2 cos 1) / (tan 2 sin 1 – tan 1 sin 2)] . (16) 

The reason this equation cannot be used is that 2 = 0° and 2 = 0° or ±180°, resulting in 

tan 2 = sin 2 = 0, causing the denominator to be zero and the value of 0 to be undefined. 

Fortunately, because the equator and the great circle meet at the points (2 = 0°,2 = 0°) and 

(0°,±180°), then 0 must be ±90°. Once again, the symmetry argument is applied as shown in 

figure 13. The figure illustrates the instance where (2 = 0°,2 = ±180°). The other point where 

the equator and the great circle meet (0°,0°) is obscured by the x-z plane. To be the point with 

maximum latitude, it must be halfway between the points where the equator and the great circle 

intersect, which is at 0 = ±90°. With two values of 0 established, it is possible to calculate 2 

values of max using equation 14, because 1 and 1 are not on the equator. It is then a matter of 

selecting the value of 0 which resulted in a positive value of max.
14

Figure 13. 2 and 2 are the intersection of the international date line and the equator. 

Alternatively, should (1,1) be at the intersection of the equator and the prime meridian or the 

equator and the international date line, equation 16 becomes undefined for the same reason it did 

when (2,2) was at the intersection of the equator and the prime meridian (or the equator and the 

international date line). And for the same reason, 0 must once again be ±90°. But this time, 

equation 11 is used to determine the two values of max. 

If 0 tends toward ±90°, then tan 0 tends toward ±infinity. This is consistent with the 

denominator of equation 16 tending toward zero. This gives us confidence that in this instance, 

our choice of 0 = ±90° was the correct one. 

3.4.3 Both (2,2) and (1,1) on the Equator 

If both (2,2) and (1,1) are on the equator, then the only great circle that can contain both 

points is the equator. In that instance, max = 0°, and 0 is undefined. The equation for the great 

circle is  = 0°, and  can be any value between –180° and +180°. Figure 14 illustrates this case. 

3.4.4 Both (2,2) and (1,1) on the Same or Opposite Meridians 

Initially, we will assume that at least one of the two points is not on the equator or at the pole. 

This means that both the numerator and denominator in equation 16 will be finite so that two 

values of 0 can be established. However, the value of max cannot be found so readily. To see 

why, we consult figure 15. 

In this case, the mathematical relationship between 1 and 2 is that they are either equal or differ 

by 180°. This causes the plane of the great circle (plane 2 in the figure) to be at a right angle to 

the equatorial plane. This means that max = 90°, causing tan max to be undefined (see equation 

14). Note that the sphere has been removed from the figure for clarity. 
15

Figure 14. Both (2,2) and (1,1) are on the equator. 

Figure 15. The two waypoints are on opposite meridians.
16

Note, too, that the value of 0 returned by equation 16 differs by ±90° from 1 and ±90° from 2, 

which, in turn, is equal to or ±180° from 1. This is consistent with the latitude 0 sitting 

midway between the latitudes of the intersection points of the equator and the great circle. In 

this case, the great circle is a meridian line, so the latitude of the intersection points of the 

meridian and the equator is equal to 1 and is ±180° from 1. However, due to the symmetry of 

the great circle and its orientation to the equator (see figure 15), it is impossible to tell which 

value of 0 corresponds to the value of max = 90°. Therefore, the equation of the great circle 

with both (2,2) and (1,1) on the same or opposite meridians will be 

  = 1, while –90° ≤  ≤ +90°, and if 1 < 0,  also = 1 + 180°: 

 if 1 > 0,  also = 1 – 180°, while –90° ≤  ≤ +90°. (17) 

In equation 17, there are only two possible values of  yet  may take on any value. Note that 

equation 17 is equally valid if (1,1) is replaced by (2,2). 

3.4.5 (2,2) or (1,1) at a Pole 

This means that 1 or 2 equals ±90°. In this case, the equation for the great circle becomes 

 = j, while –90° ≤  ≤ +90° and if j < 0,  also = j + 180°: if j > 0,  also = j – 180°, 

while –90° ≤  ≤ +90°, where j represents the number of the point that is not at the pole. Note 

that the equation is very similar to equation 17. Note too that in this case, 1 may not be equal 

to, or differ from, 2 by 180°. All that is necessary is for 1 or 2 to equal ±90°. 

3.4.6 (2,2) and (1,1) Are Antipodal 

Antipodal points are points whose coordinates are on exactly opposite sides of the sphere. This 

means that 2 =  ± 180°, and 2 = –1. Applying this to the definition for 0 as listed in 

equation 16, we find 

 0 = arc tan {[tan 1 cos (1 ±180°) – tan (–1) cos 1] / 

 [tan (–1) sin 1 – tan 1 sin (1 ± 180°)]} . (18) 

For the sine and cosine functions, changing the argument of the functions by 180° changes the 

sign of the function (see figure 9). In other words (9), 

 sin ( ± 180°) = –sin , and cos ( ± 180°) = –cos  for all . (19) 

This reduced equation 18 to 

 0 = arc tan {cos 1[tan 1 + tan (–1)] / sin 1 [tan (–1) + tan 1]} . (20)
17

The tangent function is a periodic, odd function (10) meaning that tan (–) = –tan . This causes 

equation 20 to reduce to zero divided by zero, independent of the latitude or longitude. Hence, 

we realize that when the waypoints are antipodal, they do not define a unique great circle. To 

see why, we consult figure 16. 

Figure 16. Two antipodal points do not define a unique great circle. 

We choose an arbitrary point (1,1) in plane 1, which is the radius distance R from the center 

point at the intersection of the x, y, and z axes, which is the sphere’s center point. Next, we 

choose the antipodal point (–1,1 ± 180°), also a distance R from the center point. Because the 

antipodal points are opposite each other on the sphere, a straight line can be drawn connecting 

the antipodal points and the sphere’s center point. This straight line, 2R in length, serves as the 

diameter of the circle in plane 1 in the figure. In 3-D space, however, one diameter can serve as 

the diameter for an infinite number of circles, all of which are contained in a sphere of radius R, 

as shown by the circles in planes 2 and 3 in figure 16. An infinite number of great circles can be 

rotated around their common diameter. Therefore, two antipodal points on a sphere do not 

define a great circle. Note too that antipodal points need not be on the equator. Again, the 

sphere and the equatorial plane have been omitted from the figure for clarity.
18

If two waypoints are antipodal on the Earth’s surface, the shortest distance between them will be 

on a path drawn through the poles. This is because the Earth’s shape is very slightly flattened at 

the poles due to its rotation. So on our perfect sphere model, the equation for the great circle for 

antipodal points will be the same as for two points on the same or opposite meridians discussed 

in the previous section. The equation is equation 17, which is repeated here for convenience as 

equation 21: 

  = 1, while –90° ≤  ≤ +90°, and if 1 < 0,  also = 1 + 180°: 

 if 1 > 0,  also = 1 – 180°, while –90° ≤  ≤ +90°. (21) 

Since the points are antipodal, it is equidistant to go through either the North or South Pole. For 

consistency, we chose to go through the North Pole for all antipodal waypoints. 

3.4.7 (2,2) and (1,1) at Opposite Poles 

In this case, 1 and 2 are at ±90°. The infinite number of great circles are meridians. In this 

case, we chose the equation that contains the prime meridian and the international date line. So, 

 = 0° and 180° while –90° ≤ ≤ 90°. 

3.4.8 (2,2) and (1,1) Nearly Antipodal on the Equator 

In this case, 1 = 2 = 0, and | 2 – 1| is almost, but not quite, 180°. With the waypoints not 

antipodal, yet both on the equator, the equation of the great circle connecting them can be found 

using equations derived earlier in the instance when both waypoints were on the equator. The 

great circle connecting them would be the equator. 

This would be true if the Earth were a perfect sphere. Unfortunately, it is slightly flattened at the 

poles, and when cartographic precision is needed, the Earth is modeled as an oblate spheroid 

(11), which is a special case of an ellipsoid (12), based on the WGS84 (World Geodetic System 

for 1984) model (13). In this instance, it may well be that the shortest path that connects the 

nearly antipodal waypoints passes closer to the poles than to the equator. By symmetry, the path 

could go past either the North or the South Pole. For simplicity, we decree that the path passes 

the North Pole (see figure 17).
19

Figure 17. The shortest path on the oblate spheroid connecting two equatorial waypoints 

nearly opposite each other goes closer to the poles than the equator. The plane 

provides the cross section of the oblate spheroid. 

The two waypoints can be connected with a straight line that goes through the Earth. Next, we 

superimpose a coordinate system such that the x axis is parallel to the line that connects the 

nearly antipodal waypoints. The line has the equation x = any real value, y = y0, where y0 is the 

distance between the line connecting the waypoints and the x axis, and z = 0. This puts the line 

connecting the waypoints in the x-y plane. It is also skew and perpendicular to the z axis. 

The question that arises is just how close to the pole the shortest path line passes and what the 

functional dependence of the distance is on the value of y0. To answer the question, we examine 

the cross section of the oblate spheroid in the y-z plane, viewing along the x axis in the negative 

x direction. Figure 18 shows the cross section. The cross section taken through the poles is an 

ellipse (14) with the following equation (15): 

 1 = y2

 / a2 + z2 / c2 , (22)
20

Figure 18. The y-z cross section of the oblate spheroid of figure 17. 

where a is the equatorial radius and c is the polar radius. The shortest path line on the surface of 

the oblate spheroid intersects the cross-sectional ellipse at a point fairly near the pole. The 

straight line connecting the waypoints intersects the y axis at the point (0,y0,0). These two points 

are separated by a distance D. The square of the distance D is represented by the equation 

 D2

 = (y – y0)

2

 + z2 , (23) 

which can be rendered in terms of a single variable after solving the elliptical cross-sectional 

equation (equation 20) for z2

 and substituting into equation 23 

 D2 = (y–y0)

2

 + c2 (1 – y2 / a2) . (24) 

It follows that the shortest path line along the surface of the oblate spheroid should be the line 

closest to the straight line intersecting the waypoints. Therefore, the point on the ellipse where 

the shortest mean sea level path intersects the elliptical cross section should be a point on the 

ellipse that is closest to (0,y0,0). This point should have a value for y such that the distance 

squared, D2

, is minimized. The value for y is easily found by taking the y derivative of equation 

24 and setting it to zero. The answer is 

 y = y0 / (1 – c2 / a2 ) . (25) 

The greatest value y can have is the equatorial radius a, so substituting that value in equation 25 

and solving for y0 renders the maximum value possible for y0, y0-max rendered 

 y0-max = a (1 – c2 / a2 ) . (26) 

For two points on the equator and nearly antipodal, the line connecting them on the surface of an 

oblate spheroid does not follow the equator if they are within two y0 of being 180° opposite each 

other. If they are closer than that, the closest path lies on the equator. Plugging in the WGS84 
21

values for the Earth, for a = 6,378,137.0 m and c = 6,356,752.314 245 m (13), we find 

y0-max = 42697.6727075391 m. 

To facilitate the derivation of the value of the minimum angle between two waypoints on the 

equator for which the shortest path between them is not on the equator, we define the angle max

as follows: 

 max = arc sin (y0-max / a) . (27) 

The minimum angle between two waypoints on the equator for which the shortest path is not on 

the equator then becomes 

 min = 180° – 2 max . (28) 

To see why, we consult figure 17. We suppose that we are looking in the negative z direction at 

the special case where the waypoints are located on the line y = y0-max, z = 0, and that we are 

looking from the North Pole ( = 90°) in the negative z direction and view the earth in cross 

section at the equator. We see the layout as shown in figure 19. (Note that “a” denotes the 

equatorial radius. The oblate spheroid has been removed for clarity.) 

Equation 27 describes the simple geometric relationship between y0-max, the equatorial radius a, 

and the angle max found in the figure. From the figure is evident that the angle min plus twice 

max equals a straight angle, or 180o. Solving for min produces equation 28. 

Figure 19. Cross-sectional view of the equator of figure 15 

looking in the negative z direction when y0 is at its 

maximum value and the definition of the angles 

min and max added. 
22

Applying the values of the WGS84 Geodetic System model (13) to equations 27 and 28 gives us 

the minimum separation for two waypoints on the equator that requires the path connecting them 

to go off of the equator. We conclude that if two waypoints on the equator are between 180° and 

179.232874830387° apart in longitude, the shortest path connecting them is not along the 

equator. In describing angles separating the waypoints on the equator, keep in mind the fact that 

the longitude range is between –180° and +180°, and that the functions used to derive the angle 

is periodic. Taking these facts into account, the relationship between 1 and 2, which results in 

the shortest path between them laying off of the equator, is 

 180° < 180° – |180° – |2 – 1|| < 179.232874830387° . (29) 

This is the quantitative definition of “nearly antipodal.” 

The error introduced by using a sphere to predict the shortest path between two nearly antipodal 

points on the equator can be reduced greatly by using a special process. The process involves 

using the shortest path intersection point of figure 18. By symmetry, we know that the point is 

halfway between the two waypoints on the shortest path line and midway on the straight line 

connecting the two waypoints. (This is true for any two nearly antipodal points on the equator—

the shortest path intersection point need not be in the y-z plane, as in figure 17.) The situation is 

illustrated in figure 20. 

Figure 20. Two nearly antipodal waypoints on the equator bisected by a plane that is 

perpendicular to the straight line connecting the nearly antipodal waypoints.
23

We added a plane that is defined by the origin, the North and South Poles, and the shortest path 

line intersection point (s,s). The plane is perpendicular to the straight line connecting the 

nearly antipodal points, and the plane bisects the sphere. 

Having established that s is midway between 1 and 2, we cannot simply sum 1 and 2 and 

divide by 2. The value of 1, 2, and s is limited to the interval –180° to +180°. This fact, and 

the periodic nature of the variables, requires a more careful definition of s. The complete 

definition is outlined in equation 30: 

 If | 2 – 1 | < 180°, s = ( 1 + 2 ) / 2 , 

 If | 2 – 1 | > 180°, and ( 1 + 2 ) / 2 < 0, s = [ ( 1 + 2 ) / 2 ] + 180°, 

 If | 2 – 1 | > 180°, and ( 1 + 2 ) / 2 ≥ 0, s = [ ( 1 + 2 ) / 2 ] – 180°. (30) 

We now derive the value of s, the latitude of the shortest path line intersection point. We will 

assume that the relationship between the waypoints, the straight line intersecting them, and the 

axis is as shown in figure 21. If the final mathematical expression defining s is independent of 

Cartesian variables, then the expression will be true for all longitudes. This is due to the polar 

axial symmetry of an oblate spheroid, and the fact that any ellipse produced by a cross section 

that includes both poles has the same major and minor radii. 

Figure 21. A more general view of figure 17: the straight line 

connecting the waypoints is not necessarily parallel to 

the x or y axis, yet it is perpendicular to the z axis.
24

To find the value of s, we first consider a more general view of figure 17 where the line 

connecting the waypoints is not necessarily parallel to the y or x axis, as in figure 21. Yet it is 

still perpendicular to the z axis and contained in the x-y plane. The first step in finding s is to 

find the angle . Keeping in mind the fact that the limits of 1 and 2 are between –180° and 

+180°, we see from the figure that 

 If | 2 – 1 | < 180°,  = (180°– | 2 – 1 | ) / 2 , 

 If | 2 – 1 | > 180°,  = (| 2 – 1 | –180° ) / 2 , (31) 

and 

 h = a sin  (32) 

The next step involves a more generalized view of the cross section in figure 18, as shown in 

figure 22. The value of h has already been established in equation 32. We can then find w with 

equation 25, where y0 is the quantity h, and y becomes w. 

 w = h / (1 – c2 / a2) . (33) 

Figure 22. A more general cross section of figure 18. 

In figure 22, d is the minimum distance between the straight line connecting the waypoints and 

the shortest path line intersection point, and r is the distance from the center point of the oblate 

spheroid to the ellipse. The quantity L is found from solving the equation for an ellipse using the 

known quantities w, a, and c. 

 L2

 = c2 ( 1 – w2 / a2 ) . (34)
25

Finally, the value of s is given as 

 s = arc tan (L / w) . (35) 

The final step to approximating the path between two nearly antipodal waypoints that are on the 

equator is to construct two great circles: one going from the point (0,1) to (s,s) and one going 

from point (s,s) to (0,2). The equation for the two segments takes on the form of equation 6, 

repeated here as equation 36 for convenience. 

 tan  = tan max cos ( – 0) . (36) 

Both segments will have the form of equation 36, but the values of max and 0 will be different. 

For the segment going from (0,1) to (s,s), the constants in equation 36 will be called max-1

and 0-1, while the constants corresponding to the segment from (s,s) to (0,2) will be called 

max-2 and 0-2. To find the values, we take advantage of the fact that (0,1) and (0,2) are on the 

equator. Because of this, 0-1 and 0-2 will be 90° removed from(0,1) and (0,2). So, 

 If | 2 – 1 | < 180° and 2 > 1, 0-1 = 1 + 90° 

 If | 2 – 1 | < 180° and 2 < 1, 0-1 = 1 – 90° 

 If | 2 – 1 | > 180° and 2 > 1, 0-1 = 1 – 90° and if 0-1 < –180°, 0-1 = 0-1 + 360°. 

 If | 2 – 1 | > 180° and 2 < 1, 0-1 = 1 + 90° and if 0-1 > or 180°, 0-1 = 0-1 – 360°. (37) 

Having found the value of 0-1, we can establish the value of max-1 by using equation 11. 

 max-1 = arc tan [tan s / cos (s – 0-1)] . (38) 

A similar process is used to find 0-2. 

 If | 2 – 1 | < 180° and 2 > 1, 0-2 = 2 – 90° 

 If | 2 – 1 | < 180° and 2 < 1, 0-2 = 2 + 90° 

 If | 2 – 1 | > 180° and 2 > 1, 0-2 = 2 + 90° and if 0-2 > 180°, 0-2 = 0-1 – 360°. 

 If | 2 – 1 | > 180° and 2 < 1, 0-2 = 2 – 90° and if 0-2 < or 180°, 0-2 = 0-2 + 360°. (39) 

As with max-1, max-2 is found using equation 11. 

 max-2 = arc tan [tan s / cos (s – 0-2)] . (40)

Note that max-1 and max-2 must always be kept in the range 0° to 90°.
26

The shortest line has now been approximated as two distinct great circle segments. One has the 

equation tan  = tan max-1 cos ( – 0-1), for the segment with  between 0-1 and s, while the 

other has the equation tan  = tan max-2 cos ( – 0-2), for the segment with  between 0-2 and s. 

This may appear to be a crude approximation for the actual path connecting two nearly antipodal 

waypoints that lie on the equator. However, a previous investigation (1, pp 17–19, 27) 

discovered that navigation software provided by the U.S. Army Topographic Engineering Center 

called FORWARD and INVERSE (16) either could not calculate any path between two nearly 

antipodal waypoints on a WGS84 ellipsoid, or it predicted a path along the equator 

(1, pp 17–19, 27). Although this module uses a perfect sphere as its basis, it will be able to 

produce a path closer to the truly shortest path between nearly antipodal points. 

3.4.9 Summary of Waypoint Placement Problems 

Table 1 summarizes the difficulty of waypoint placement discussed in the previous sections. The 

difficulty arises when the waypoints are used to calculate the parameters max and 0, which 

define the great circle. Methods to circumvent their difficulties (discussed in the previous 

sections) are also summarized in the table. 

Table 1. Summary of waypoint placement problems and their solutions. 

Waypoint Placement Problem Solution

Nearly antipodal on equator. Path connecting them is 

nonequatorial; sphere model predicts 

equatorial path. 

Calculate intermediate point with 

equations 30–35; use equations 

36–40 to define two great circle 

connecting segments. 

Waypoints at opposite poles. Infinite number of meridians all 

equally valid. 

Chose to use the prime meridian. 

Waypoints antipodal. Antipodal points do not define great 

circle. 

Use great circle that passes through 

North Pole. 

Waypoints on same or opposite 

meridian. 

Equation 11 is undefined. Use great circle that contains 

meridians. 

One waypoint at the pole. Equation 10 is undefined. Use great circle that contains 

meridians. 

Both waypoints on the equator. Equation 10 is undefined. Use the equator as great circle. 

One waypoint at intersection of 

equator and prime meridian or 

international date line. 

Equation 10 is undefined. Assign 0 = ±90°, use equation 11 

or 14 to calculate max, use value of 

0, which gives max > 0. 

3.5 Value of the Earth’s Average Radius 

In a previous work (1, pp 17–19, 27) it was established that modeling the Earth as a perfect 

sphere provides a basis for interwaypoint distance calculation. The distances calculated would 

be accurate enough to serve the needs of the software. So, the question arises as to how to 

establish the value of the radius of the perfect sphere.
27

Our approach is to pick a radius such that if a latitude and longitude were randomly chosen on a 

WGS84 oblate spheroid, sea level would have a 50% chance of being higher than the average 

radius and a 50% chance of being lower than the average radius. This entails adding every 

possible radius and dividing by the number of radii. Since there are an infinite number of radii in 

an oblate spheroid, we integrate all possible radii over all space, then divide by the integration 

over all space to define the average radius Rave as where r(,) is the radius of the latitude  and 

longitude , and dS is the solid angle differential, also a function of latitude and longitude. 

 

. (41)

The solid angle differential in spherical coordinates is written dS = cos d d (17). (Note that 

the reference’s differential is written as dS = sin  d d. This is because the reference’s angle 

is measured down from the North Pole. In our coordinate system, the reference angle  is 

measured up from the equator.) To obtain the radius of an oblate spheroid, we consider the 

oblate spheroid’s general equation in Cartesian coordinates (18), which is 1 = (x2 + y2) / a2 + z2 / 

c

2

, where a is the equatorial radius and c is the polar radius. Rewriting this in spherical 

coordinates using equations 1, 2, and 3, and using r (the dummy variable in lieu of R), we have 

 1 = ( r2

 cos 

2  cos 2 + r2

 sin 2  cos 

2  ) / a2

 + r2 sin 2  / c2 . (42) 

When we apply the trigonometric identity cos2

 + sin2 = 1, solve for r, and change the 

integration limits to radians from degrees, using equation 42, equation 41 becomes 

 . (43)

Evaluating the denominator is straightforward: the value is 4. Evaluating the variable  in the 

numerator is also straightforward: the value is 2. This is because an oblate spheroid is 
28

symmetric about the z axis, making it independent of . With those variables evaluated, 

equation 43 then becomes 

 

. (44) 

To make the integral wieldy, a variable substitution is made such that u = sin . Taking the 

differential du = cos  d, changing the limits, and doing some algebraic manipulation renders 

the average radius Rave as 

 

. (45)

Evaluating this integral (19, 20), applying the limits of integration, and performing some 

algebraic manipulation render the average radius as 

 (46)

because the equatorial radius a is greater than the polar radius c. After applying the WGS84 

values for the equatorial radius a and the polar radius c as 6,378,137.0 m, and c = 

6,356,752.314245 m, respectively (13), the value of Rave is 6370994.40182752 m, which varies 

by <0.0008% from the value offered by the Chemical Rubber Handbook (21), which was not 

derived using the WGS84 oblate spheroid model. 

4. Results 

The formulae derived in spherical coordinates from the analysis of section 3 is next applied to 

obtaining the results of algorithms and formulae needed to describe platform motion. After a 

presentation of terms used to describe moving platforms, the equations of motion are derived 

assuming constant acceleration and constant deceleration up to and down from the platform’s 

cruising velocity. Next, an algorithm is presented for adding, removing, and modifying 

waypoints. Finally, assuming that the snapshot times are known, the location of the platform 

between waypoints at the snapshot times is obtained. 
29

4.1 Definition of Moving Platform Variables 

For all moving platforms, we assume that the velocity vs. time profile takes a consistent structure 

between two waypoints. The relevant variables and constants used to describe the profile are 

displayed in figure 23. As a platform approaches the jth waypoint (labeled “Waypoint j:” in the 

figure), it is traveling at a constant maximum cruise velocity, vcru. This maximum cruise velocity 

is a function of the kind of platform being modeled (high-mobility multipurpose wheeled vehicle 

[HMMWV] vs. tank vs. unmanned aerial vehicle [UAV] vs. dismounted Soldier, etc.) and the 

kind of terrain the platform is crossing (paved road vs. unpaved road vs. smooth terrain vs. rocky 

terrain, etc.) assuming that it is a ground platform. The value for the vcru is read from a database 

of platform performance specifications provided by Army development centers. If the platform 

is part of a squad or convoy, vcru can be assigned the value of the maximum cruising velocity of 

the slowest vehicle in the squad or convoy to ensure that the platforms remain in formation. 

Figure 23. Variables that describe moving platforms as they travel between waypoints. 

Note that vcru is not the maximum velocity of which the platform is capable; rather, it is a typical 

speed the platform obtains when it is in travel mode. For example, for a dismounted Soldier vcru

might be the brisk walking pace of 3 mph, even though a Soldier doing double time can travel at 

nearly twice that rate. We assign the value that is typical of a dismounted Soldier in travel mode.
30

When the platform approaches the jth waypoint, it begins to slow down with an acceleration of 

–a0. The minus sign ensures that a0 is assigned a positive value. This value, like vcru, is also a 

function of the kind of platform as well as the terrain the platform crosses. It, too, will be read 

from a database with information supplied by Army development centers. And as with the 

variable vcru, this will not be the maximum value of which the platform is capable but rather a 

typical value for the platform. 

When the platform arrives at the jth waypoint, it will do so with a velocity of vj at the arrival 

time ta,j. The waypoint location is specified in latitude (j), longitude (j), and elevation (Ej). 

The elevation is the distance above the average sea level for a perfectly spherical Earth, which is 

the average radius Rave from the Earth’s center point. The radius (notated as “R” in figure 2) for 

the jth waypoint is Rave + Ej. The latitude and longitude are specified by the user, while the 

elevation Ej is the elevation of the ground on which the platform is traveling (from a database 

with information provided by the Digital Terrain Elevation Data [DTED] [22]) plus the antenna 

height above the ground. If the platform is not a ground platform, the elevation Ej is specified by 

the user. 

When the platform departs the jth waypoint for the j+1th waypoint, it will begin to accelerate 

with the value +a0 at the departure time td,j. This is the same value with which it decelerated 

before reaching the waypoint, only in the opposite direction of the deceleration –a0. (The orange 

hash in the figure that crosses the blue path is the point where the platform begins accelerating or 

decelerating.) So the value for a0 must be chosen to represent the platform’s normal acceleration 

and deceleration. The departure velocity and the arrival velocity at the jth waypoint must be the 

same, vj. If the platform does not stop at the waypoint but simply passes through it, then it 

leaves the waypoint with the same speed with which it arrived. If it does stop, then it must leave 

the waypoint from a stop, so the arrival and departure speeds are zero. If the platform does not 

stop at the waypoint, then the arrival time ta,j and departure time td,j must also be the same. If the 

platform does stop, then it is possible for the platform to loiter at the waypoint for a waiting, or 

loiter, time tl,j = td,j – ta,j. 

The platform accelerates until it reaches the cruise velocity vcru. It then stops accelerating, 

traveling at a constant speed of vcru. When it approaches the j+1th waypoint, it begins to 

decelerate at –a0, until it reaches the j+1th waypoint with speed vj+1 at time tj+1,a. As with the jth 

waypoint, if the platform has a speed of zero when it reaches the j+1th waypoint, it may loiter for 

time tl,j+1 = td,j+1 – ta,j+1 before it leaves the waypoint at time td,j+1 when it accelerates at a0 to reach 

vcru on its way to the j+2th waypoint. 

The path the platform travels between two waypoints is the most expeditious path between two 

points on a perfect sphere: an arc of a great circle (23). As the platform passes the waypoint, it 

changes direction instantaneously so as to follow the great circle path to the next waypoint. 

Because the direction changes, the platform’s velocity changes instantly. Yet, in an attempt to 
31

maintain a greater degree of realism in modeling platform motion, the platform’s speed is treated 

as constant as it passes through a waypoint. 

4.1.1 Mathematical Description of Platform Motion 

Figure 24 shows graphs of acceleration, speed, and distance vs. time of platform motion between 

waypoints. To keep the motion as simple as possible, we define the acceleration graph as a 

simple step function: a0 at the departure time for waypoint j td,j. At time t1, it becomes 0, and at 

t2 it becomes –a0, until it arrives at waypoint j+1 at time ta,j+1. At all waypoints, the acceleration 

is 0. 

Figure 24. Motion profile for moving platform between waypoints. 

Integrating the acceleration graph yields the speed vs. time graph. The speed begins at vj at 

waypoint j, then increases at a constant acceleration rate of a0 until time t1 when the acceleration 

goes to zero, and the speed reaches vcru. The speed remains at vcru until time t2, when it decreases 

at a deceleration rate of –a0 to become vj+1 at time ta,j+1. 

To find the distance covered, we integrate the speed graph in time. At time td,j the platform 

moves away from the waypoint x(j,j,Ej), such that the distance increases with time and the rate 

of distance change increases in time until time t1. From time t1 to time t2, the distance increases 

at a constant rate in time. From time t2 until the platform reaches waypoint j+1 at time ta,j+1, the 

distance is still increasing in time but at a decreasing rate.
32

4.1.2 Mathematical Relationship Between the Motion-Defining Variables 

To begin, we define the acceleration of the platform in terms of the time variable t according to 

the graph in figure 24. 

 At t = td,j, a(t) = 0 . (47) 

 At td,j < t < t1, a(t) = a0 . (48) 

 At t = t1, a(t) = a0 . (49) 

 At t1 < t < t2, a(t) = 0 . (50) 

 At t = t2, a(t) = –a0 . (51) 

 At t2 < t < ta,j+1, a(t) = –a0 . (52) 

 At t = ta,j+1, a(t) = 0 . (53) 

To find the speed dependence in time, we integrate equations 48, 50, and 52 in time. The 

constant of integration will be determined by making the speed at times td,j and ta,j+1 equal to vj

and vj+1, respectively. We will also use the fact that the speed from figure 24 shows that from 

(and including) times t1 to t2, the speed is vcru. We will also assume that the speed is continuous 

for all time. From the boundary condition we know that 

 at t = td,j v(t) = vj . (54) 

Integrating equation 48 in time, we know that 

 at td,j < t < t1, v(t) = a0 t + C , (55) 

where C is the constant of integration. Applying the boundary condition of equation 54, we 

evaluate C, and with some algebra we find that 

 at td,j < t < t1, v(t) = a0 ( t –td,j) + vj . (56) 

Evaluating the expression at the first inflection point t1, we see that the speed should be vcru from 

figure 24, so 

 at t = t1, v(t) = a0 ( t1 – td,j) + vj = vcru . (57) 

Between the inflection times t1 and t2, the speed remains constant at vcru, so that 

 at t1 < t < t2, v(t) = vcru , (58)
33

and it follows from figure 24 that 

 at t = t2, v(t) = vcru , (59) 

giving us another boundary condition. Integrating equation 47 in the time window t2 < t < ta,j+1, 

we find that 

 at t2 < t < ta,j+1, v(t) = –a0 t + C , (60) 

where applying equation 59 to evaluate the integration constant C and rearranging the variable 

shows us that 

 at t2 < t < ta,j+1, v(t) = vcru – a0 ( t – t2 ) . (61) 

To complete the speed profile in time between the waypoints, we use the speed at the j+1th 

waypoint to find that 

 at t = ta,j+1, v(t) = vcru – a0 ( ta,j+1 – t2 ) = vj+1 . (62) 

To obtain the distance D(t) traveled by the platform from waypoint j as a function of time, 

equations 56, 58, and 61 are integrated in time. The result for equation 56’s integration is 

 at td,j < t < t1, D(t) = a0 t2

 / 2 – a0 td,j t + vj t + C . (63) 

The value of the integration constant C is found by knowing that D(t) = 0 at time t = td,j. 

Solving for C and rearranging the variables, we find 

 at td,j < t < t1, D(t) = a0 ( t2 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( t – td,j ) . (64) 

The distance the platform obtains from the waypoint j at the inflection time t1 is found by 

applying t1 to the variable t in equation 64, where 

 at t = t1, D(t) = a0 ( t1

2

 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( t1 – td,j ) . (65) 

Distance D(t) must increase from the value outlined in equation 65. It will increase with the 

constant speed vcru as shown in equation 58. So, the integration of equation 58 is added to 

equation 65 so that in the next time segment, where t1 < t < t2, the distance formula is given 

 at t1 < t < t2, D(t) = a0 ( t1

2

 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( t1 – td,j ) + vcru t + C . (66) 

To find the integration constant C, we use the fact that at time t1, equation 66 equals equation 64. 

The result is that 

 at t1 < t < t2, D(t) = a0 ( t1

2

 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( t1 – td,j ) + vcru (t – t1 ) . (67) 

The distance of the platform D(t) from waypoint j at the inflection time t2 is 

 at t = t2, D(t) = a0 ( t1

2

 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( t1 – td,j ) + vcru (t2 – t1 ) . (68)
34

The distance will increase from the point where the platform is located at time t2. The rate at 

which the distance D(t) will increase in time is outlined in equation 61, so that adding the 

integration of equation 61 shows that 

 at t2 < t < ta,j+1, D(t) = a0 ( t1

2

 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( t1 – td,j ) + vcru (t2 – t1 ) + 

vcru t – a0 ( t2

 – 2 t2 t ) / 2 + C . (69) 

Evaluation of the integration constant C is accomplished by equating equation 69 with equation 

68 for time t = t2, so the result is that 

 at t2 < t < ta,j+1, D(t) = a0 ( t1

2

 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( t1 – td,j ) + vcru (t2 – t1 ) + 

vcru ( t – t2 ) – a0 ( t2 – t2

2

 ) / 2 + a0 t2 ( t – t2 ) . (70) 

At time t = ta,j+1, equation 70 gives the distance between waypoints j and j+1. For this reason, we 

notate the distance between them as Dj, evaluate equation 70 making t = ta,j+1, and gather some 

terms together to simplify the equation. After rearranging some of the variables, we find that 

 at t = ta,j+1, Dj = a0 ( t1 – td,j )2 / 2 + vj ( t1 – td,j ) + vcru (ta,j+1 – t1 ) – a0 (ta,j+1 – t2 )2 / 2 . (71) 

Bringing back equation 57 and equation 62, and after applying some algebra to isolate and 

number the time variables, equations 71–73 are three independent equations showing the 

relationship between nine variables. 

 t1 – td,j= ( vcru – vj ) / a0 . (72) 

 ta,j+1 – t2 = ( vcru – vj+1 ) / a0 . (73) 

By combining equations 71–73, it is possible to eliminate the inflection times t1 and t2, which 

delineate the times when the platform will change speed. Doing so produces an equation that 

defines the distance between the jth and j+1th waypoint Dj in terms of the departure time of the 

jth waypoint td,j, the arrival time at the j+1th waypoint ta,j+1, the platform speed at the jth 

waypoint vj, the platform speed at the j+1th waypoint vj+1, the platform cruising speed between 

waypoints vcru, and the platform acceleration/deceleration rate a0. This is equation 74. 

 Dj = vcru (ta,j+1 – td,j) – (vcru – vj+1)

2

 / (2 a0) – (vcru – vj)

2

 / (2 a0) . (74) 

4.1.3 Restrictions on Input Data – General 

When the user specifies the location of the j = 1st waypoint, the latitude and longitude will need 

to be input. The latitude can range from –180° to +180°, while the longitude ranges from –90° to 

+90°. If the waypoint is for a ground platform, the elevation is read from software containing 

elevation data, like the DTED database (22). If the platform is not a ground platform, the user 

inputs the elevation data, restricted by the ground elevation, and the platform’s service ceiling, 

which is read from a database of platform specifications. By default, the corresponding arrival 
35

time ta,1 will be set to 0 s, since this is the beginning of the simulation. Also by default, the 

platform speed will be set to 0 m/s if it is a ground platform or the minimum air speed if it is an 

air platform. The minimum air speed will be read from a database with specification data for the 

platform. The user will then set the loiter time, defining the departure time for the first waypoint 

td,1, which, for j = 1, will be equal to the loiter time. The cruise speed vcru and the acceleration a0

are read from the platform specification database. 

Next, the user inputs the latitude and longitude for the second waypoint. The elevation is also 

input if the platform is an air platform, otherwise the elevation is read from a DTED-like 

database. This defines the distance between the waypoints D1. Consulting equation 74 for j = 1, 

all variables but two have been defined: ta,2 and v2. If the user defines the time ta,2, then 

waypoint speed v2 is also defined. If the user defines the waypoint speed v2, then by equation 

74, ta,2 is also defined. If the user specifies a loiter time at waypoint two (i.e., making it so that 

ta,2 is not equal to td,2), then v2 must be equal to 0 m/s. The result in that case is that ta,2 will be 

calculated according to equation 74, independent of any user desire.

In the general case (i.e., for j >1), the distance Dj is defined when the user inputs the waypoint 

latitudes, longitudes, and elevations (or the elevations are read from an elevation database). vcru

is also read from a platform specification database, as is a0. tj and vj were defined when the 

calculations for the last waypoints were done, leaving the user free to choose either ta,j+1, which 

defines vj+1 by equation 74, or vj+1, which defines ta,j+1. Should the user decide that the platform 

loiter at the j+1th waypoint, then vj+1 must be zero, defining the arrival time ta,j+1 and giving the 

user no flexibility insofar as having variables to define. 

Equation 74 is rewritten to show the arrival time and the speed at the j+1th waypoint as a 

function of the other variables. Equation 73 was obtained by simply solving equation 75 for 

ta,j+1. 

 ta,j+1 = td,j + [ 2 a0 Dj + (vcru – vj+1)

2

 + (vcru – vj)

2

 ] / (2 a0 vcru ) . (75) 

Equation 76 was derived by solving equation 74 for the speed at the j+1th waypoint. When the 

variable vj+1 was isolated in equation 74, it became evident that there were two possible solutions 

for vj+1, as it was necessary to take a square root to find vj+1. The two solutions found were vcru

plus an expression under the radical and vcru minus an expression under the radical. We take the 

minus case because the radical’s value is always positive, since a platform’s speed cannot exceed 

vcru. Once the speed at the j+1th waypoint is calculated with equation 76, 

, (76)

and the j+1th waypoint arrival time ta,j+1 is established, equations 72 and 73 can be used to find 

the inflection times t1 and t2. 

When the last waypoint is assigned, we assume that the speed of the platform becomes zero at 

the final waypoint. Hence, the user has the freedom to assign the location of the final waypoint, 
36

while the arrival time is determined by equation 75 for the final point. Should the simulation 

stop before the platform reaches the final waypoint, the final time is taken as an interpolation 

point, and the final latitude, longitude, and elevation are calculated, treating the final time like an 

interpolation time. Should the simulation end time occur after the time the platform reaches the 

final waypoint, the platform will simply treat the time between the arrival time at the final 

waypoint ta,j+1 and the simulation end time as a loiter time at the waypoint. 

4.1.4 Distance Calculation 

Using the spherical model, we next calculate the distance between two points notated 

(Rave+Ej,j,j) and (Rave+Ej+1,j+1,j+1). Rave is the average radius of the Earth; Ej, j, and j are 

the jth’s waypoint elevation, latitude, and longitude; and Ej+1, j+1, and j+1 are the j+1th’s 

waypoint elevation, latitude, and longitude. The points are shown with the great circle they 

define, and a wedge containing them and the center point in figure 25. The equator is in blue, the 

great circle in red, and the wedge formed by the waypoints and the Earth’s center point is shown 

in green. Note that the waypoints are not on the great circle but are elevated above it by Ej and 

Ej+1, respectively. 

Figure 25. Two waypoints and the Earth’s center point 

form a wedge. 
37

To find the distance between the waypoints, we must remove the sphere and take a vantage point 

normal to the wedge formed by the two waypoints and the Earth’s center. This is shown in 

figure 26. Rave is the Earth’s average sea level radius, Ej and Ej+1 are the jth and j+1th waypoint 

elevations, j and j+1 are the latitudes of the waypoints, and j and j+1 are the longitudes.  is 

the angle formed at the Earth’s center point between the lines connecting the waypoints with the 

center point. Dj is the distance along the platform’s path between the jth and j+1th waypoint, 

while  is the line-of-sight distance between the waypoint projections on the sea level surface of 

the Earth. 

Figure 26. The wedge from figure 25 with the Earth removed for clarity. 

First, we find the relationship between the angle  and the line-of-sight distance . From the 

two back-to-back right triangles, we see that 

  / Rave = 2 sin (  / 2 ) . (77) 

Written another way, equation 77 becomes 

  = 2 arc sin [  / ( 2 Rave ) ] . (78) 

Rave is known but  is not. To find  we employ the well-known distance formula in a 3-D, 

rectilinear coordinate system (24). 

 . (79) 

Substituting equations 1, 2, and 3, into 79, then employing equation 78, the term Rave cancels to 

reveal  in terms of the latitudes and longitudes of the two waypoints. 

  = 2 arc sin { [( 1 / 2 ) [(cos j+1 cos j+1 – cosj cos j)

2

 + 

(cos j+1 sin j+1 – cos j sin j)

2

 + (sin j+1 – sin j)

2

]

1/2} . (80)
38

Note that the sign of  is positive because we are taking the positive root of equation 78. 

Next, to find the arc length Dj we employ the definition of arc length in planer motion (25), 

which is appropriate to figure 26 because the wedge is contained in a single plane. 

 

. (81) 

We rewrite equation 81 so that the distance becomes S, a continuous variable as a function of the 

continuous variable . We keep the variable R to be the distance from the Earth’s center point to 

the platform. So, for any point on the path S the platform traverses, the angle it makes with the 

previous waypoint  as measured at the Earth’s center point, is 

 

. (82)

We know that the elevation at a waypoint Ej is related to the distance from the Earth’s center 

point Rj in that 

 Rj = Rave + Ej . (83) 

We now use the variables E and R to represent the continuous elevation and continuous distance 

of the platform from the Earth’s center point on the path between waypoints. Updating equation 

83, we see that 

 R = Rave + E . (84) 

Differentiating equation 84 in time, it becomes self-evident that the rate of change of elevation 

dE/dt is equal to the rate of change of radius dR/dt, as stated in equation 83. 

 dR / dt = dE / dt . (85) 

Now we establish the relationship between the change in time of path length S and the change in 

time of the distance of the platform from the Earth’s center point R. For an aerial platform like a 

fixed wing aircraft or a helicopter, we assume that the rate of climb or descent is directly 

proportional to its airspeed along the path S. For a ground platform, the equivalent assumption is 

that the elevation change is proportional to the distance traversed. Put another way, we assume 

that the slope of the ground between waypoints is constant. Expressed quantitatively, in both the 

aerial platform and the ground platform case, 
39

 dE / dt = dR / dt = k dS / dt , (86) 

where k is a constant whose value will be determined in due course. 

Next, we wish to establish the relationship between the platform’s distance from the Earth’s 

center R, the path length between waypoints S, and the angle the platform makes at the Earth’s 

center between the platform’s location and the previous or jth waypoint. Multiplying equation 

86 by the time differential dt and integrating, we see that 

 R = k S + C , (87) 

where C is the constant of integration. To evaluate it, we know that when S = 0, the platform is 

at the jth waypoint, so the distance from the Earth’s center to the waypoint will be Rj = Rave + Ej. 

This allows us to evaluate the integration constant in equation 87. After the constant’s 

evaluation, equation 87 is rendered 

 R = k S + Rj = k S + Rave + Ej . (88) 

To evaluate k, it will be necessary to differentiate equation 88 by . 

 dR / d = k dS / d. (89) 

Note that the differential of the constant terms is zero. Next, we define the antideferential F() 

such that dF() / d = [ ( dR / d )2 + R2 ]

1/2. Rewriting equation 82 in terms of F(), we arrive at 

the distance traveled from the last waypoint to the platform. 

 S = F() – F(0) . (90) 

Differentiating equation 90 and substituting for the term dS / d from equation 89 (keeping in 

mind that F(0) is a constant, so its derivative is zero), we find 

 dR / d = k [ ( dR / d )2

 + R2 ]

1/2. (91) 

Squaring equation 91 and grouping the differentials together, then taking the square root, we 

arrive at a conventional first-order differential equation. 

 (dR / R) = ± d k / (1 – k2

)

1/2 . (92) 

As the platform moves, the value of  increases. If the elevation increases, then the value of k 

must be positive in equation 92. If the elevation decreases, then the value of k must be negative 

in equation 92, since we assume that the quantity (1 – k2

)

1/2 is always positive. Consider the 

constant k as listed in equation 88. If the value of the elevation of the platform increases with 

distance, k is positive. Likewise, as it decreases with distance, k is negative. So, we will simply 

let the sign of k change to indicate the elevation change of the platform: if k is positive, the 

platform is ascending. If k is negative, the platform is descending. To maintain consistency in 

equation 92, we take the positive root. 
40

Solving equation 92 and taking the natural antilogarithm, the equation becomes 

 R = C exp [ k / (1 – k2

)

1/2] , (93) 

where R is the distance of the platform to the Earth’s center point,  is the angular distance along 

the path, exp [ ] signifies the value e (roughly 2.718281828459 [26]) raised to the power of the 

expression inside the brackets, and C is the constant of integration. To find C, we use initial 

conditions: when  = 0, the right half of equation 93 becomes C, which is equal to the left half, 

which is R = R1 = Rave + E1. 

To find the value of k, we evaluate equation 93 when the platform reaches the j+1th point and 

equation 80 to find the angular distance between the waypoints. 

 Rj+1 = Rave + E j+1= Rj exp [ k / (1 – k2

)

1/2] = (Rave + Ej ) exp [ k / (1 – k2

)

1/2] . (94) 

Taking the natural logarithm of equation 94 and solving for k, we find that 

 k = ln [ (Rave + Ej+1) / (Rave + Ej) ] / { 2 + { ln [(Rave + Ej+1) / (Rave + Ej)]}2 }1/2 , (95) 

transforming k to a known value. Next, we incorporate equation 88 and use the values of the 

variables when the platform reaches the j+1th waypoint. 

 Rj+1 = Rave + Ej+1 = k Dj + Rj = k Dj + Rave + Ej . (96) 

After solving for the distance the platform travels between the waypoints Dj and eliminating the 

average radius of the Earth, we find the distance in terms of k and the elevations at each 

waypoint. 

 Dj = (Ej+1 – Ej) / k . (97) 

Combining equations 94 and 96 and eliminating Rj+1, the radius at the j+1th waypoint, and 

algebraically rearranging some variables, it is possible to find an alternative expression for the 

distance the platform travels. 

 Dj = (Rj / k ) {exp [ k / (1 – k2)

1/2] – 1} = 

[ (Rave + Ej )/ k ] {exp [ k / (1 – k2)

1/2] – 1} . (98) 

The constant k denotes the change of elevation divided by the distance traveled. So, the possible 

ranges of k can be expressed as –1 ≤ k ≤ 1. Equation 98 becomes undefined if k = ±1, but 

equation 97 will still produce a finite value for Dj. However, both equations 97 and 98 are 

undefined if k = 0. It is possible to take the limit as k approaches zero in equation 98 and use the 

result to compute Dj. We find 

 Dj = Rj  = (Rave + Ej)  (99)
41

consistent with the formula for the length of a circle’s arc (27)—the radius times the angle in 

radians. Since k = 0, the platform does not change altitude, making the path the arc of a circle. 

Having established the platform path distance Dj, it is now possible to calculate the arrival time 

ta,j+1 using equation 75 or the arrival velocity vj+1 using equation 76. 

4.1.5 Restrictions on Input Data – Maximum Speed, Minimum Time 

The maximum speed a platform can attain between waypoints is vcru, the cruising speed. To 

arrive at the j+1th waypoint in the minimum amount of time, the speed profile would have to be 

as shown in figure 27. The profile differs from the one displayed in figure 24 only in that the 

final speed remains at the cruising velocity, thus minimizing the arrival time ta,j+1. Additionally, 

the two inflection times t1 and t2 have been replaced with a single inflection time ti. It follows 

that the arrival speed vj+1 = vcru. Since figure 27 is a special case of figure 24, the equation 

describing the arrival time ta,j+1 is a special case of equation 75. Equation 75 is restated here as 

equation 100. 

 ta,j+1 = td,j + [ 2 a0 Dj + (vcru – vj+1)

2

 + (vcru – vj)

2

 ] / (2 a0 vcru ) . (100) 

To find the minimum arrival time, vcru is set equal to vj+1 in equation 100, and the arrival time is 

notated as tmin. 

 tmin = td,j + [ 2 a0 Dj + (vcru – vj)

2

 ] / (2 a0 vcru ) . (101) 

Figure 27. The motion profile for maximum arrival speed and minimum arrival time.

For some speeds at waypoint j and some distances between waypoints j and j+1, the platform 

may not have time to accelerate to the speed vcru. In that case, the motion profile for speed vs. 

time looks like figure 28. The distance covered Dj is the average speed (vj + vmax)/2 times the 

time to travel between waypoints tmin – td,j, or 

 Dj = ( tmin – td,j ) (vj + vmax) / 2 . (102)
42

Figure 28. The maximum speed if the distance does not allow enough time for Vmax to equal Vcru.

The relationship between the time traveled between waypoints tmin – td,j, the acceleration a0, and 

the speeds vj and vmax is 

 ( vmax – vj ) = a0 ( tmin – td,j ) . (103) 

After eliminating the variables tmin and td,j using equations 102 and 103, then solving for vmax, 

we find the maximum speed in this case is 

 vmax = ( 2 a0 Dj + vj

2

)

1/2. (104) 

Using equation 104 to eliminate vmax from equation 103, and solving for tmin, we find the 

corresponding minimum arrival time at the j+1th waypoint is 

 tmin = td,j + [(2 a0 Dj + vj

2

)

1/2

 – vj ] / a0 . (105) 

So, the maximum speed at the j+1th waypoint vj+1 the user can select for the platform is vmax, as 

defined by equation 104. If the value of equation 104 does not exceed vcru, then the 

corresponding minimum arrival time ta,j+1 the user may select is tmin, defined by equation 105. If 

the value for vmax as defined by equation 104 exceeds the platform cruising speed vcru, then the 

maximum arrival speed the user can select for vj+1 is vcru with a minimum arrival time ta,j+1, for 

the platform is defined by equation 105. 

4.1.6 Restrictions on Input Data – Minimum Speed, Maximum Time 

The minimum speed with which a platform can arrive at the j+1th waypoint is generally (but not 

always) zero. This is a special case of the velocity profile of figure 24, shown here as figure 29. 

Equations 75 and 76 govern the platform’s motion (with equations 72 and 73 defining the 

inflection times t1 and t2) with vj+1 set to zero. The corresponding maximum time is 

 tmax = td,j + [ 2 a0 Dj + vcru2 + (vcru – vj)

2

 ] / (2 a0 vcru ) (106) 

and the inflection times are 

 t1 = td,j + ( vcru – vj ) / a0 (107) 

and 

 t2 = tmax – a0 / vcru . (108)
43

Figure 29. The motion profile for zero arrival speed and maximum arrival time. 

4.1.7 Restrictions on Input Data – Peak Speed Less Than vcru. 

For some speeds at waypoint j and some distances between waypoints j and j+1, the platform 

may not be able to achieve the cruising speed vcru. In that case, the speed profile resembles 

figure 30. To find the value of vpk, we invoke equations 71–73. But in this case, there is only 

one inflection point ti, and at no point between the jth and j+1th waypoints does the platform 

travel at a constant speed. The arrival speed vj+1 is zero. With these considerations in mind, 

equation 71 for the distance can be rewritten as 

 Dj = (vpk + vj) (ti – td,j) / 2 + vpk (ta,j+1 – ti) / 2 . (109) 

Figure 30. The speed of the platform if it never reaches vcru. 

Equations 72 and 73 relate the inflection times and the start and ending times with the 

acceleration and speed. They are 

 ti – td,j = ( vpk – vj ) / a0 (110) 

and 

 ta,j+1 – ti = –vpk / –a0 . (111) 

We use equations 110 and 111 to eliminate the departure time td,j, the inflection time ti, and the 

arrival time ta,j+1 from equation 109, rendering it as 

 Dj = (vpk

2

 – vj

2

) / (2 a0) + vpk

2

 / (2 a0) . (112) 
44

Upon solving for vpk, we find 

 vpk = [ (2 a0 Dj + vj

2

) / 2 ]1/2 . (113) 

To find the corresponding maximum time tmax, we eliminate the inflection time ti in equations 

110 and 111. 

 td,j + ( vpk – vj ) / a0 = ta,j+1 – vpk / a0 . (114) 

The final step is to substitute equation 115 for vpk into equation 116 and solve for ta,j+1 and define 

it as tmax. 

 tmax = td,j + { 2 [ (2 a0 Dj + vj

2

 ) / 2 ]1/2 – vj } / a0 . (115) 

4.1.8 Restrictions on Input Data – Insufficient Time for Platform to Reach Zero 

For some speeds at waypoint j and some distances between waypoints j and j+1, the platform 

may not have time to come to a full stop. In that case, the motion profile for speed vs. time looks 

like figure 31. In this case, the arrival velocity is the minimum velocity calculated with equation 

116 but with the sign for the acceleration a0 changed. 

 vmin = (vj

2

 – 2 a0 Dj)

1/2. (116) 

Figure 31. The minimum speed if the platform does not have time to reach zero speed at the arrival waypoint. 

Likewise, the expression for the maximum arrival time resembles equation 105, but the sign for 

the acceleration is changed. 

 tmax = td,j + [vj – (vj

2

 – 2 a0 Dj )1/2] / a0 . (117) 

So, the minimum speed at the j+1th waypoint vj+1 the user can select for the platform is zero, but 

only if the quantity vj

2

 – 2 a0 Dj (see equation 109) is negative. The corresponding maximum 

time tmax is defined by equation 106, with 107 and 108 giving the values for the inflection times 

t1 and t2. If the quantity vj

2

 – 2 a0 Dj is positive, the value for the minimum arrival time is set by 

equation 109. The corresponding maximum time tmax is calculated with equation 117. 
45

Should the user wish to have a loiter time, the arrival speed va,j+1 must be set to zero by default. 

However, should the user desire a loiter time but the condition that the quantity vj

2

 – 2 a0 Dj is 

>0, then the arrival speed is set according to equation 117, and the loiter time is automatically set 

to zero. 

4.1.9 Inputting Speed or Time Data 

To facilitate computing the data needed for interpolation, we first consider the speed profile as 

shown in figure 32. It defines the top speed vtop of the platform between waypoints for the 

profile shown in the figure, which happens at the inflection time tt. The figure is a special case 

of the speed profile of figure 24. Likewise, the quantitative description of the figure is a special 

case of equations 72, 73, 75, and 76. To make use of those equations, the expression for vtop

must be derived. 

Figure 32. The quantity vtop and tt are defined graphically. 

Suppose the user specifies the arrival velocity vj+1. The distance between the waypoints will be 

 Dj = (vtop + vj) (ti – td,j) / 2 + (vj+1 + vtop) (ta,j+1 – ti) / 2 . (118) 

The inflection time ti and start and ending times td,j and ta,j+1 are related to the acceleration and 

velocities vj, vj+1 by equations 119 and 120. 

 tt – td,j = ( vtop – vj ) / a0 . (119) 

 ta,j+1 – tt = ( vj+1 – vtop ) / –a0 . (120) 

After substituting equations 119 and 120 into equation 118 to eliminate all the time variables, we 

find 

 Dj = (vtop

2

 – vj

2

) / ( 2 a0 ) – (vj+1

2

 – vtop

2

) / ( 2 a0 ) . (121) 

Solving for vtop renders 

 vtop = [ (2 a0 Dj + vj+1

2

 + vj

2

 ) / 2 ]1/2 . (122)
46

If the user specifies the arrival time instead, we first use equations 119 and 120 to eliminate tt

and solve for vj+1. 

 v j+1 = 2 vtop – vj – a0 ( ta,j+1 – td,j ) . (123) 

We next solve equation 115 for vj+1. 

 v j+1 = (2 vpk

2

 – vj

2

 – 2 a0 Dj )1/2 . (124) 

Finally, we equate equations 123 and 124 to eliminate vj+1. After solving for vtop, we see that 

 vtop = vj + a0 (ta,j+1 – td,j ) – [ vj a0 (ta,j+1 – td,j ) – a0 Dj + ( 1 / 2 ) a0

2

 ( ta,j+1 – td,j )2 ]1/2. (125) 

When solving to find equation 125, the negative case was selected on the basis of calculating test 

values for the variables. 

Lastly, we restate equations 72, 73, 75, and 76 with vtop substituted for vcru, and isolate the 

inflection times t1 and t2. 

 t1 = td,j + ( vtop – vj ) / a0 . (126) 

 t2 = ta,j+1 – ( vtop – vj+1 ) / a0 . (127) 

 ta,j+1 = td,j + [ 2 a0 Dj + (vtop – vj+1)

2

 + (vtop – vj)

2

 ] / (2 a0 vtop ) . (128) 

 vj+1 = vtop – {2 a0 [ vtop (ta,j+1 – td,j) – Dj ] – (vtop – vj)

2

 }1/2 . (129) 

If the user selects the velocity vj+1, we calculate the top speed vtop using equation 122. If vtop is 

less than vcru, vtop retains its value. If vtop is greater than vcru, then vtop is set equal to vcru. We 

then use equation 128 to calculate ta,j+1, and equations 126 and 127 to calculate the inflection 

point times t1 and t2. 

If the user decides to select the j+1th arrival time ta,j+1, we calculate the top speed vtop using 

equation 125. If vtop is less than vcru, vtop retains its value. If vtop is greater than vcru, then vtop is 

set equal to vcru. We then use equation 129 to calculate vj+1, and equations 126 and 127 to 

calculate the inflection point times t1 and t2. 

4.1.10 Summary of Inputting Waypoint Data 

Figure 33 shows the overall process of inputting waypoint data. The process shown in the figure 

is inclusive of equations 47 to 129, although it does not include the restrictions outlined in table 

1. The rectangular boxes represent processes to be executed. The diamond boxes are questions 

with yes or no answers. The answers determine the calculations, tasks, and the order in which 

the tasks and calculations are performed. The green letters indicate the “Start” box, and the 

variables defining the platform’s motion are in blue. Violet shows the indices, amber indicates 

equation numbers, and red shows when it is necessary to “Do Interpolation.” 

At the beginning, the user identifies the platform type and the number of waypoints (the 

recording index), while the software retrieves the acceleration rate a0 and the cruising speed vcru 
47

Figure 33. Overview of the waypoint data input process.
48

based on the platform type. If the platform is not a ground platform, the user inputs latitude, 

longitude, and elevation. In the case that the platform is a ground platform, the software reads 

the elevation data from the DTED database (22). The first point motion data, v1 and ta,1, are set 

to zero. The user inputs the loiter time, and the index j is set to 1. Implicit in these processes is 

the assumption that every platform has at least two distinct waypoints. The case where a 

platform has only one waypoint will be treated later. 

If the index indicates that the next waypoint is the last, a special process begins. As with the first 

waypoint, the user inputs the latitude, longitude, and (if the platform is not a ground platform) 

elevation. When this data is used, the distance from the jth waypoint is calculated. If the distance 

is not sufficient to allow the platform to stop (vj

2

 is not >2 a0 Dj [remember, this is the final 

waypoint and the platform must arrive here with a speed of zero]) while it decelerates at a0, the 

deployment module assigns the minimum possible speed at the waypoint vj+1. It then calculates 

the corresponding arrival time ta,j+1, departure time td,j+1 (remember, there is no loiter time because 

the speed is not zero), and inflection time ti. Lastly, it forces the user to input another waypoint (in 

this case, t1 = t2, so the i represents the subscripts 1 and 2 since there is only one inflection point). 

This process continues until the platform traverses enough distance to come to a halt. If there is 

enough distance for the platform to come to a halt, then the speed vj+1 is set to zero and the greatest 

speed between the waypoints vpk is calculated. If vpk is greater than the specification vcru, then vpk

is set to the value of vcru; t1, t2, and ta,j+1 are calculated; and the interpolated points are calculated. If 

vpk is less than the specification vcru, then vpk remains unchanged, and ti (again, the i subscript 

represents 1 or 2 since in this case t1 = t2) is calculated with ta,j+1, and the interpolated points are 

also calculated. 

If the next waypoint isn’t the final point, the user is asked to input the latitude, longitude, and (in 

the event the platform is not a ground platform) the elevation. The distance between the 

previous waypoint and the current waypoint Dj is calculated. If there is not enough distance for 

the platform to reach the cruising velocity vcru, then the maximum velocity vmax is calculated; 

otherwise, vmax is set equal to vcru. The minimum corresponding arrival time tmin is calculated. If 

there is not enough distance for the platform to reach zero velocity, the minimum velocity vmin is 

calculated, and the loiter time is set to zero and the corresponding value of tmax calculated. Then, 

vmin, tmax, vmax, and tmin are displayed as guidance. Otherwise, vmin is set equal to zero, with the 

values tmin, tmax, vmin, and vmax displayed as guidance. If vmin is zero, the user is prompted to 

input the loiter time. 

If the loiter time is zero, then final data for the platform motion at the waypoint is inputted, either 

as an arrival time ta,j+1 or waypoint speed vj+1. The deployment module calculates the maximum 

possible speed between the current and previous waypoint vtop and compares it to the cruise 

speed vcru. If vtop is greater than vcru, then vtop is set equal to vcru. When vtop is used, the 

inflection points t1 and t2 are calculated, and the departure time td,j+1, after which the index j is 

increased by 1, and the index are tested to see if the next waypoint is the final waypoint.
49

4.2 Modifying Input 

Next, we assume that a string of waypoints has been inputted and recorded by the deployment 

module, but that the user wishes to modify them. Four processes are possible: moving a 

waypoint, adding a waypoint, removing a waypoint, and changing either the arrival time ta,j+1, 

waypoint velocity vj, or if the waypoint velocity vj is zero, the loiter time. 

4.2.1 Moving a Waypoint 

Figure 34 shows the sequence for moving a waypoint. It is similar to, but somewhat more 

abbreviated than, the process for inputting a waypoint. After starting the process, the user 

identifies the platform and waypoint to be moved. If this is a ground platform, the user inputs 

the new latitude and longitude while the software provides the elevation data from the DTED 

database (22). If it is not a ground platform, the user provides the latitude, longitude, and 

elevation. If this is the first waypoint, the index is immediately changed because the waypoint 

speed and the initial time are zero and must not be changed. The loiter time can be changed by 

another process. Implicit in all these processes is the assumption that every platform has at least 

two distinct waypoints. The case where a platform has only one waypoint will be treated later. 

If this is the final waypoint for the platform, the distance Dj is calculated, and if there is enough 

distance for the platform to stop (2 a0 Dj > vj

2

), the final speed vj+1 is assigned to zero and the 

peak speed vpk is calculated. If the platform peak speed vpk exceeds the cruising speed vcru, then 

vpk is set to vcru, and the arrival times ta,j+1, t1, and t2 (see figure 29) are calculated and the 

interpolation redone. If the platform is not greater than vcru, then the single inflection time ti is 

calculated (in this case, t1 = t2, so the i represents the subscripts 1 and 2, since there is only one 

inflection point) along with the arrival time ta,j+1 (see figure 30), after which the interpolation is 

redone. Should the distance Dj not be sufficient to allow the platform to stop, vmin is calculated 

and the speed vj+1 is assigned the value vmin; the arrival time ta,j+1 is calculated, along with the 

inflection time ti (in this case, t1 = t2, so the i represents the subscripts 1 and 2); tloiter is assigned 

to 0; and the index j is increased by one, as is the recording index (which is the number of 

waypoints for this platform). Then data is input for the next waypoint (latitude, longitude, and 

elevation if it is not a ground platform, and latitude and longitude with elevation read from a 

database if it is a ground platform), and the distance Dj calculated. The distance is then tested to 

see if the platform has enough room to stop. If it is not, the user is prompted for another point 

until the platform has distance enough to stop. 

If this is not the final waypoint, the distance Dj is calculated. If there is enough distance Dj for 

the speed to reach vcru, then the maximum speed of the platform between waypoints vmax is 

assigned the value vcru. If there is not enough distance, then the maximum speed the distance 

allows (equation 104) is assigned to vmax. Next, the speed for this waypoint assigned previously 

(vj+1) is tested to see if it is greater than the vmax. If it is, then vj+1 is assigned the value vmax and 

if it is not, vj+1 remains unchanged. Then, the distance is tested to see if the platform has enough 

space to reach zero velocity. If it does, then vmin is assigned the value of 0. If it does not, then 
50

Figure 34. Overview of the process of moving a waypoint. 
51

the minimum possible velocity vmin is assigned the lowest speed the platform can reach (vj

2

 – 2 a0

Dj)

1/2 by equation 116. vj+1 is next tested to see if it is less than vmin. If it is, it is assigned the 

value of vmin. If it is not, it remains unchanged. Next, vtop is calculated by equation 122. If vtop

is greater than the cruise speed vcru, vtop is assigned the value of vcru. If it is not, it remains 

unchanged. Next, if the waypoint velocity vj+1 is nonzero, the loiter time is set to zero. If vj+1 is 

zero, then the loiter time remains unchanged. Finally, ta,j+1, t1, t2, and td,j+1 are calculated, and the 

index j is increased to check the next waypoint. 

4.2.2 Adding a Waypoint 

The flowchart for adding a waypoint is shown in figure 35. Portions of the figure resemble the 

process of moving a waypoint, shown in figure 34. After identifying the platform type and new 

index number, the user inputs the longitude, latitude, and elevation for an aerial platform, or just 

the latitude and longitude for a ground platform. The indices for the waypoints that follow the 

new waypoint in time are incremented by one, as is the recording index that records the number 

of waypoints. If the added waypoint is intended to be the first waypoint, the software assigns the 

waypoint speed and arrival time as zero. The user is prompted for the loiter time, the departure 

time td,1 is set to the loiter time, and the index is increased by one. The next waypoint is tested to 

see if it is the final waypoint. 

If the added waypoint is not intended to be the first waypoint, it is tested to see if it is the final 

waypoint. Should this be the final waypoint, the distance Dj is calculated. If there is enough 

distance for the platform to stop (2 a0 Dj > vj

2

), the final waypoint velocity vj+1 is assigned to zero 

and the peak speed vpk is calculated. If the peak speed vpk has time to get to vcru, then vpk is 

assigned the value of vcru, and the arrival times ta,j+1, t1, and t2 are calculated (see figure 29) and 

the interpolation redone. Otherwise (see figure 30), vpk remains unchanged and ti is calculated 

(the i subscript represents 1 or 2, since in this case t1 = t2) with ta,j+1 and the interpolated points 

are recalculated. If there is not enough distance for the platform to stop, vmin is calculated and 

the speed vj+1 is assigned the minimum possible value for speed vmin; tmax is calculated and 

assigned to the arrival time ta,j+1; tloiter is assigned to 0; and the index j is increased by one, as is 

the recording index (which is the number of waypoints). The user is prompted to input the 

latitude, longitude, and elevation for an aerial platform, or only latitude and longitude for a 

ground platform. The distance Dj to the new platform is calculated, and the cycle of new 

platforms is inputted until it is able to stop. 

If this is not the final waypoint, then guidance must be calculated to let the user select the new 

waypoint arrival time ta,j+1 or waypoint velocity vj+1. This is first accomplished by calculating 

the distance to the new waypoint Dj. If there is enough distance for the platform to reach vcru, 

(2 a0 Dj + vj

2

)

1/2, then vmax is assigned the value of vcru, and the corresponding tmin is calculated. 

If there isn’t enough distance, then maximum speed vmax is calculated, with the corresponding 

tmin.
52

Figure 35. Overview of the process of adding a waypoint. 
53

Next, the guidance values for vmin and tmax are calculated. If there is enough distance Dj for the 

platform to stop (2 a0 Dj – vj

2

 > 0), then vmin is assigned a value of zero and the corresponding 

value of tmax is calculated. The formula to use for finding tmax depends on whether there is 

enough distance for the speed to reach vcru. If there is [ vcru < (2 a0 Dj + vj

2

 )1/2 ], then tmax is 

calculated using equation 106. If there is not, tmax is calculated with equation 115, followed by 

the display of vmin, tmax, vmax, and tmin to the user as guidance, who inputs the loiter time. If there 

is not enough distance Dj for the platform to stop (2 a0 Dj – vj

2

 < 0), then vmin is calculated with 

equation 116, and the corresponding value of tmax is calculated with equation 117 and displayed 

to the user as guidance while the loiter time is set to zero. 

If the loiter time is not zero, the waypoint speed vj+1 is set to zero. The highest speed possible for 

the platform between waypoints vtop is calculated with equation 122, and if it exceeds vcru, vtop is 

set to vcru, and ta,j+1 is computed with equation 128. If the loiter time is zero, the user is given the 

option of entering the arrival time ta,j+1 or the waypoint velocity vj+1. If the user chooses to enter 

ta,j+1, then vtop is calculated using equation 125, and if it exceeds vcru, vtop is set to vcru, and vj+1 is 

computed using equation 129. If the user chooses to enter vj+1, then vtop is calculated with 

equation 122, and if it exceeds vcru, it is set to vcru, and ta,j+1 is computed with equation 128. Then 

t1 and t2 are calculated with equations 126 and 127, td,j+1 is calculated, and the waypoint index is 

increased by 1. 

If the index is for the final waypoint, then it is treated as the final waypoint as outlined 

previously. If it is not, it and all the remaining waypoints will have to be checked for selfconsistency. To do this, the distance Dj is calculated, and from this, vmax is calculated with 

equation 104. If vmax exceeds vcru, vmax is set to vcru, otherwise it retains the same value. The 

waypoint speed, vj+1, is compared to vmax; if it is greater than vmax, it is set to vmax. If it is not, it 

retains its value. Similarly, if 2 a0 Dj – vj

2

 is positive, vmin is set equal to zero, or the square root 

of its negative if 2 a0 Dj – vj

2

 is negative. Then vj+1 is compared to vmin and set equal to vmin if 

vj+1 is less than vmin. A new variable, vtop is computed with equation 122. If vtop is greater than 

the cruising speed, vcru, vtop is set equal to vcru. If vj+1 is greater than zero, the loiter time is set to 

zero. Then ta,j+1, t1, t2, and td,j+1 are calculated, and the waypoint index is increased by 1. 

4.2.3 Removing a Waypoint 

Removing a waypoint is a shorter procedure and is illustrated in figure 36. The user first 

identifies the index and type of platform for which the waypoint must be removed. After reading 

the acceleration factor ao and the top cruising speed vcru from a database, the remaining indices 

for the waypoints for this platform have 1 subtracted from them, and the recording index has 1 

subtracted from it. 

If the waypoint removed was the first, the former second waypoint becomes the new first 

waypoint. Then, the waypoint velocity and arrival time are set to zero. The user is prompted for 

the loiter time, and the index is increased by one. The next waypoint is tested to see if it is the 

last waypoint. 
54

Figure 36. Overview of the process of removing a waypoint. 
55

If this is the final waypoint, the distance Dj is calculated, and if there is enough distance for the 

platform to stop (2 a0 Dj > vj

2

), the final waypoint velocity vj+1 is assigned to zero and the peak 

speed vpk is calculated with equation 113. If vpk is greater than vcru, then vpk is assigned the value 

of vcru, the arrival times ta,j+1, t1, and t2 are calculated (see figure 29), and the interpolation is 

redone. Otherwise (see figure 30), vpk remains unchanged and ti (the i subscript represents 1 or 

2, since in this case t1 = t2) and ta,j+1 are calculated and the interpolated points recalculated. If 

there is not enough distance for the platform to stop, the speed vj+1 is assigned the minimum 

possible value for speed vmin after calculating vmin and the arrival time ta,j+1, while tloiter is 

assigned to 0 and the index j is increased by one, as is the recording index (which is the number 

of waypoints). The user is prompted to input the latitude, longitude, and elevation for an aerial 

platform, or only latitude and longitude for a ground platform. The distance Dj to the new 

platform is calculated and the cycle of new platforms inputted until it is able to come to a stop. 

If it is not the final waypoint, it and all the remaining waypoints will have to be checked for selfconsistency. To do this, the distance Dj is calculated, and from this, a test is done to see if the 

maximum speed (2 a0 Dj + vj

2

)

1/2 exceeds vcru. If it does, vmax is set to vcru. The waypoint speed 

vj+1 is compared to vmax; if it is greater than vmax, it is set to vmax. Similarly, if 2 a0 Dj – vj

2

 is 

positive, vmin is set equal to 0, or if 2 a0 Dj – vj

2

 is negative, vmin is set equal to the square root of 

vj

2 – 2 a0 Dj. Then vj+1 is compared to vmin and set equal to vmin if vj+1 is less than vmin. A new 

quantity, vtop, is calculated, and if it exceeds the cruising speed vcru, it is set to vcru. If vj+1 is nonzero, the loiter time is set to zero. Then ta,j+1, t1, t2 and td,j+1 are calculated, and the waypoint 

index is increased by 1. 

4.2.4 Changing Waypoint Data 

The final process involves changing the waypoint speed, the arrival time ta,j+1, or the loiter time 

of a single waypoint. Figure 37 summarizes this. The user first identifies the index and type of 

platform for which the waypoint data is to be changed. After reading the acceleration factor ao

and the top cruising speed vcru from a database, the waypoint index is tested to see if it is the 

initial point. If it is, the process prompts the user to input a new loiter time. The arrival and 

departure times are then changed for the remaining waypoints by the difference of the old and 

new loiter times in the first point. It is unnecessary to change the other velocities because the 

first point velocity must remain zero. The process then stops. Should this not be the initial 

waypoint, the index is tested to see if it is the final waypoint. If it is, the process stops since the 

last waypoint must have a waypoint velocity vj+1 of zero, and it is unnecessary to change the 

arrival time. 

If the index is not for the final waypoint, the distance from the last waypoint Dj is recalculated. 

If there is enough distance for the platform to reach vcru, (2 a0 Dj + vj

2

)

1/2, then vmax is assigned 

the value of vcru, and the corresponding tmin is calculated with equation 101. If there isn’t enough 

distance, then maximum speed vmax is calculated with equation 104, with the corresponding tmin

calculated with equation 105. Next, if there is not enough distance for the platform to reach zero
56

Figure 37. Overview of the process of changing the waypoint data. 
57

velocity, the minimum velocity vmin is calculated with equation 116, with the corresponding tmax

calculated with equation 117, and the loiter time is set to zero, and tmin, tmax, vmin, and vmax are 

displayed to the user as guidance. Otherwise, vmin is set equal to zero, and the maximum 

corresponding arrival time tmax is calculated with either equation 106 (if the maximum speed vmax

is vcru) or with equation 115 (if the maximum speed vmax is vcru). The values tmin, tmax, vmin, and 

vmax are displayed to the user as guidance. If vmin is zero, the user is prompted to input the loiter 

time. 

The loiter time is then tested. If the loiter time is not zero, the waypoint speed vj+1 is set to zero. 

The highest speed possible for the platform between waypoints vtop is calculated with equation 

122, and if it exceeds vcru, it is set to vcru, and ta,j+1 is computed with equation 128. If the loiter 

time is zero, the user is given the option of entering the arrival time ta,j+1 or the waypoint velocity, 

vj+1. If the user chooses to enter ta,j+1, then vtop is calculated with equation 125, and if it exceeds 

vcru, it is set to vcru, and vj+1 is computed with equation 129. If the user chooses to enter vj+1, then 

again vtop is calculated with equation 122, and if it exceeds vcru, it is set to vcru, and ta,j+1 is 

computed with equation 128. Then t1, t2, and td,j+1 are calculated using equations 126 and 127, 

and the waypoint index is increased by 1. 

The index is again tested to see if it corresponds to the final waypoint. If this is the final 

waypoint, the distance Dj is calculated, and if there is enough distance for the platform to stop (2 

a0 Dj > vj

2

), the final waypoint velocity vj+1 is assigned to zero and the peak speed vpk is 

calculated with equation 113. If the peak speed vpk exceeds vcru, then vpk is assigned the value of 

vcru, and the arrival times ta,j+1, t1, and t2 are calculated with equations 75, 107, and 108 (see 

figure 29) and the interpolation is redone. Otherwise (see figure 30), vpk remains unchanged and 

ti (the i subscript represents 1 or 2, since in this case t1 = t2) and ta,j+1 are calculated with 

equations 110 and 111, and the interpolated points are recalculated. If there is not enough 

distance for the platform to stop, the speed vmin is calculated with equation 116 and the speed vj+1

is assigned the value for vmin, the arrival time ta,j+1 (which equals tmax) is calculated with equation 

117, tloiter is assigned to 0, and the index j is increased by one, as is the recording index (which is 

the number of waypoints). The user is prompted to input the latitude, longitude, and elevation 

for an aerial platform, or only latitude and longitude for a ground platform. The distance Dj to 

the new platform is calculated and the cycle of new platforms inputted until the platform is able 

to come to a stop. 

If it is not the final waypoint, it and all the remaining waypoints will have to be checked for selfconsistency. To do this, the distance Dj is calculated, and from this, a test is done to see if the 

maximum speed (2 a0 Dj + vj

2

)

1/2 can exceed vcru. If it can, vmax is set to vcru; otherwise, vmax is 

calculated with equation 104. The waypoint speed, vj+1, is compared to vmax—if it is greater than 

vmax, it is set to vmax. Otherwise vj+1 remains unchanged. Similarly, if 2 a0 Dj – vj

2

 is positive, 

vmin is set equal to 0, or (vj

2

 – 2 a0 Dj)

1/2 if 2 a0 Dj – vj

2

 is negative, and the loiter time is set to 

zero. Then vj+1 is compared to vmin and set equal to vmin if vj+1 is less than vmin. Otherwise, vj+1

remains unchanged. A new variable, vtop is calculated with equation 122. If it exceeds vcru, it is 
58

set to vcru. Otherwise, it remains unchanged. If vj+1 is nonzero, the loiter time is set to zero. Then 

ta,j+1, t1, t2, and td,j+1 are calculated with equations 126–128, the waypoint index is increased by 1, 

and the departure time td,j+1 is calculated. 

4.2.5 If a Platform Has Only One Waypoint 

We now consider the input process and the four modification processes (moving a waypoint, 

adding a waypoint, removing a waypoint, and reassigning values to the waypoint) when a 

platform has only one waypoint assigned to it. When a platform has only one waypoint, the 

waypoint velocity (v1) is assigned to zero. The arrival time ta,1 is also assigned to zero. The 

departure time td,1 is equal to the loiter time, which equals the mission time. If a platform has 

only one waypoint, we assume that it does not move but remains stationary for the entire run. As 

a result, the interpolation process for this platform need only copy the initial waypoint latitude, 

longitude, and elevation for every snapshot time to the snapshot database. 

Inputting data for a single waypoint platform is relatively straightforward and is set out in figure 

38. The user identifies the type of platform that requires a single waypoint, and the default index 

(as well as the recording index) is set to 1. The arrival time ta,1 and the waypoint velocity v1 are 

set to zero. The loiter time becomes the mission time, so that the departure time td,1 is set equal 

to the loiter time. If the platform is a ground platform, the user inputs the latitude and longitude, 

while the elevation is read from a DTED database. If it is not a ground platform, the user also 

inputs the elevation. The interpolation is done by copying the latitude, longitude, and elevation 

database to the snapshot database for this platform. 

Figure 38. The process of inputting data for a platform with one waypoint. 

Moving a waypoint (when a platform has a single point) only requires a reassignment of a new 

latitude and longitude. Once the user identifies the platform, the index is assumed to be 1. If it 

is a ground platform, the user inputs the new latitude and longitude. If it is not, the user inputs 

the new latitude, longitude, and elevation. Since that is all that is required, the process ends. 

This brief process is outlined in figure 39. 
59

Figure 39. The process of inputting data for changing the location of a single waypoint platform. 

Adding a point is a bit more involved. After identifying the platform on which to add a 

waypoint, the user identifies whether the added point comes before the previously defined point 

in time (the new point has an index of j = 1) or after it (the new point has an index of j = 2). The 

platform acceleration a0 and cruising speed vcru are read from a database, and the recording index 

is set to 2. If it is a ground platform, the user inputs the new latitude and longitude. If it is not, 

the user inputs the new latitude, longitude, and elevation. The two points then (if necessary) 

have their indices modified so that the first point in time has an index of 1 and the second 

waypoint in time has an index of 2. Since we now have two distinct waypoints, the process 

continues as outlined in figure 35. The entire process is laid out in figure 40. The verbiage 

describing the remaining part of the process is identical to the verbiage used to describe figure 

35. For the sake of brevity, it will not be repeated here. 

Removing a data point for a platform with only one data point is tantamount to deleting the 

platform. So instead of outlining a separate process to remove the only waypoint to a platform, 

the process will be to simply delete the platform. Likewise, reassigning waypoint data for a 

platform with only one waypoint is unnecessary since there is only one waypoint. The single 

waypoint must have the waypoint velocity v1 and arrival time ta,1 equal to zero, so the data can 

not be changed. Consequently, this process will not be included in the deployment module. 

4.3 The Interpolation Process 

Once the latitudes, longitudes, and elevations have been established for the set of waypoints for a 

platform (with their corresponding arrival and departure times either assigned directly by the 

user or calculated from assigned waypoint velocities and loiter times), then it is possible to 

calculate the platform’s latitude, longitude, and elevation for a selected snapshot time from the 

interpolation time ti. The process is to first calculate the interpolated distance Di from the most 

recent waypoint the platform encountered to the platform’s position at ti. The interpolated 

distance Di makes it possible to obtain the interpolated elevation Ei. The interpolated angle i

that the interpolated location of the platform at time ti makes with the location of the last 

waypoint at the Earth’s center then makes it possible to calculate the interpolated longitude i

and latitude i after making a coordinate transformation to simplify the use to the quantity i. 

Each of these processes will be explained in detail.
60

Figure 40. The process of adding a waypoint to a platform that has one 

waypoint. 
61

4.3.1 Calculate the Interpolated Distance 

Assuming that we know the interpolated time ti, the first task is to see if any of the waypoint 

loiter periods (the time between the arrival time tj,a and departure time td,j) bound the interpolated 

time ti. Put another way, is ta,j < ti < td,j? If it is, then the interpolated latitude, longitude, and 

elevation for the platform are set equal to the latitude, longitude, and elevation of the waypoint. 

This ends the interpolation process for the interpolated time ti. (Note that the last waypoint loiter 

period is the span between the final waypoint arrival time, tN,a and the ending time of the run.) 

If the interpolated time ti does not occur during the waypoint loiter time, then it occurs during the 

transit time between points. The interpolation process then becomes more involved. To 

calculate the interpolated distance Di, it is first necessary to compare it with the times that mark 

the inflection points t1 and t2 of the platform profile of figure 24. As a review, we recount how 

the inflection points are calculated. Step 1 is to calculate vtop from equation 125, restated here as 

equation 130.

 vtop = vj + a0 (ta,j+1 – td,j ) – [ vj a0 (ta,j+1 – td,j ) – a0 Dj + ( 1 / 2 ) a0

2

 ( ta,j+1 – td,j )2 ]1/2 . (130) 

If vtop exceeds vcru, then vtop is set equal to vcru. Otherwise, equation 130 defines vcru. Equations 

126 and 127 (repeated here as equations 131 and 132) are then used to establish the inflection 

times t1 and t2. 

 t1 = td,j + ( vtop – vj ) / a0 . (131) 

 t2 = ta,j+1 – ( vtop – vj+1 ) / a0 . (132) 

To calculate the distance Di from the most recently departed waypoint, it is necessary to 

established the relationship between the departure time td,j, the inflection times t1 and t2, the 

arrival time at the next waypoint ta,j+1, and the interpolation time ti. The next three equations 

(133–135) are taken from equations 64, 67, and 70. 

 If td,j < ti < t1, Di = a0 ( ti

2

 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( ti – td,j ) . (133) 

 If t1 < ti < t2, Di = a0 ( t1

2

 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( t1 – td,j ) + vcru (ti – t1 ) . (134) 

 If t2 < ti < ta,j+1, Di = a0 ( t1

2

 – td,j

2

 )/ 2 – ( a0 td,j – vj ) ( t1 – td,j ) + vcru (t2 – t1 ) + 

vcru ( t i – t2 ) – a0 ( ti

 2 – t22

 ) / 2 + a0 t2 ( ti – t2) . (135) 

Now that we have the interpolated distance Di, we can calculate the interpolated elevation, Ei. 

4.3.2 Calculate the Interpolated Elevation 

If the interpolated time ti does not occur during the waypoint loiter time, it is possible to 

calculate the interpolated elevation Ei as a result of one of two circumstances. First, the two 

waypoints bounding the platform in time have the same elevation. In that case, Ei = Ej = Ej+1.
62

But should Ej not equal Ej+1, the interpolated elevation must be found using a two-step process. 

First reconsider equation 97 by eliminating then solving for the elevation constant k: 

 k = (Ej+1 – Ej) / Dj . (136) 

Note that this calculation is made with the distance between the waypoints Dj, not the 

interpolated distance Di. To obtain the interpolated elevation Ei, we use equation 97 again but 

substitute the interpolated distance Di (calculated in the previous section) for the interwaypoint 

distance Dj, and the interpolated elevation Ei for the arrival waypoint elevation Ej+1. Rearranging 

the terms we find 

 k Di + Ej = Ei . (137) 

With the interpolated distance Di and the interpolated elevation Ei, we may next find the 

interpolated angle between the line connecting the most recently encountered waypoint and the 

center of the Earth, and the line connecting the platform’s position at the interpolated time ti and 

the center of the Earth. 

4.3.3 Calculate the Interpolated Angle i

Again, we assume that the interpolated time ti does not occur during the waypoint loiter time. 

And again, there are two circumstances we need to consider in the calculation of the interpolated 

angle i: the elevations of the waypoints bounding the interpolated point are either equal or 

different. 

First, if the elevations of the waypoints bounding the interpolated point in time are equal, we 

make use of equation 99, substituting the interpolated angle i for the angle the waypoints make 

at the Earth’s center , and the distance to the interpolated point Di for the distance between the 

waypoint Dj. When we solve for i 

 Di / (Rave + Ej) = i . (138) 

Second, if the elevations of the waypoints bounding the interpolate point in time are not equal, 

we use equation 94 instead. Making use of the interpolated elevation Ei and the constant k, and 

substituting the interpolated elevation Ei for the waypoint elevation Ej+1 and solving for i, we 

find 

 i = [(1 – k2

)

1/2 / k] Ln [(Rave + Ei) / (Rave + Ej )] . (139) 

Ln denotes the natural log function. Note that if the platform is descending, the elevation 

constant k will be negative (see equation 136). If the platform is descending, then Ej will be 

greater than Ei, causing the natural log of the expression to be negative. We assume that we take 

the positive root for [(1 – k2)

1/2, resulting in a positive interpolation angle i. Should the 

platform be ascending, k and the natural log will be positive, again resulting in a positive 

interpolation angle i. Note that when the elevation does not change (the situation described by 

equation 138), the interpolation angle i is positive because the distance interpolated Di, the 
63

Earth’s average radius Rave, and the waypoint elevation Ej are always positive. The fact that the 

interpolation angle i is always positive will play a vital role in determining the interpolated 

longitude and latitude.

4.3.4 Calculate the Interpolated Latitude and Longitude 

To calculate the interpolated longitude and latitude involves doing a coordination transformation 

such that the two waypoints bounding the interpolated point in time will be on the transformed 

equator. Hence, the interpolated point will also be on the transformed equator, which will be 

calculated by using the interpolated angle i. By reversing the process, we can obtain the 

interpolated longitude and latitude. 

The transformation will be done as a two-part process. The first part involves keeping the 

latitude constant but changing the longitude. This process is illustrated in figure 41. The new x 

axis (denoted x’ in the figure) begins at the Earth’s center and intersects the point on the equator 

where the line of constant latitude 0 crosses the equator. The y’ axis begins at the Earth’s 

center, is in the equatorial plane, and is 90° from the x’ axis. The axis is occulted by the Earth in 

the figure, although the original y axis (also in the equatorial plane and beginning at the Earth’s 

center) is not occulted and is plainly shown in the figure. The z and z’ axis are the same and start 

at the Earth’s center and progress through the North Pole. 

Recall that the parameters max and 0 describe the equation of a great circle and are derived from 

the latitudes and longitudes of the two bounding waypoints in equations 10 and 11, restated here 

as equations 140 and 141. 

 0 = arc tan [(tan j cos j+1 – tan j+1 cos j) / (tan j+1 sin j – tan j sin j+1)] . (140) 

 max = arc tan [tan j+1 / cos (j+1 – 0)] = arc tan [tan j / cos (j – 0)] . (141) 

Recall, too, that this assumes that the argument for arc tan is finite in equation 140, and that we 

take the positive value for max for the value of 0 calculated in equation 140. Additionally, the 

point (max,0) in the nonprimed coordinate system is the point of highest latitude on the great 

circle. Since we want the line of constant latitude that goes through this point to be the prime 

meridian in the primed system, the transformed coordinate will be (max’,0°). 

So, the effect of the first part of the transformation will be simply to change the longitude by the 

quantity 0 while allowing the latitude to remain unchanged. A way to visualize this 

transformation is to remove the Earth and look along the z axis in the –z direction. Figure 42 

illustrates the transformation. From the figure, it is easy to obtain the transformation 

mathematically for not only the interpolated points, but for any point that includes the waypoints 

bounding the interpolated point: 
64

Figure 41. The first step in transforming the interpolated point to a new coordinate system. 

Figure 42. The first step in transforming the 

interpolated point to a new 

coordinate system with the Earth 

removed viewing in the –z 

direction. 
65

 j’ = j . (142) 

 j+1’ = j+1 . (143) 

 i’ = i . (144) 

 j’ = j – 0 . (145) 

 j+1’ = j+1 – 0 . (146) 

 i’ = i – 0 . (147) 

We must make sure that the three transformed longitudes j’, j+1’, and i’ are in the range –180° 

to +180°. If any of them are less than –180°, we add 360°. If any of them exceed +180°, we 

subtract 360°. 

The second part of the transformation is more complicated. It involves defining a set of doubleprimed coordinates such that the primed and double-primed y axis are the same. Furthermore, 

the x” axis begins at the center of the sphere and goes through the point (max,0) in the original 

system and (max,0) in the single-primed system. This line intersects the double-primed equator 

and the prime meridian. The transformation is illustrated in figure 43. 

Figure 43. The second step in transforming the interpolated point to a new coordinate system: the y’ 

axis acts as the axel so that y’ = y”.
66

To obtain the quantitative relationship between the primed and double-primed coordinate 

systems, we view along the y’/y” axis in the positive y’/y” direction with the Earth removed for 

clarity. Figure 44 shows this situation. The fact that the waypoints and the interpolated points 

lie on the equator in the double-primed coordinate system simplifies the transformation. The 

next few equations will pertain to any point on the double-primed equator. As the figure shows, 

any point on the equator has an x” coordinate and a y” coordinate. Because the point is not the 

transformed double-primed equator, it has no z” coordinate. Note that y’ = y”, so that 

relationship is trivial. Concentrating on the relationship between the x’’ coordinate and the z’ 

and x’ coordinate, we see from the figure that 

 z’ = x” sin max (148) 

and 

 x’ = x” cos max . (149) 

Figure 44. The second step in transforming the interpolated 

point to a new coordinate system with the Earth 

removed viewing in the +y’ direction. 

Making use of equations 1–3 for the relationship between spherical and rectangular coordinates, 

and using the Rave, the Earth’s average radius, and the elevation E, we may rewrite equations 148 

and 149 for any point on the great circle as 

 ( Rave + E ) sin ’ = ( Rave + E ) cos ” cos ” sin max (150) 

and 

 ( Rave + E ) cos ’ cos ’ = ( Rave + E ) cos ” cos ” cos max , (151)
67

and the fact that the points on the double-primed equator have” = 0 and we may divide by the 

common factor ( Rave + E ) renders the relationships between the primed and double-primed 

latitudes and longitudes as 

 sin ’ = cos ” sin max (152) 

and 

 cos ’ cos ’ = cos ” cos max . (153) 

Substituting the relationships between the original and the double-primed coordinates as outlined 

in equations 142 to 147, and solving for the double-primed longitude ”, we find for the jth 

waypoint, the j+1th waypoint, and the interpolated point such that 

 j” = arc cos ( sin j / sin max ) . (154) 

 j” = arc cos [ cos j cos ( j – 0) / cos max ] . (155) 

 j+1” = arc cos ( sin j+1 / sin max ) . (156) 

 j+1” = arc cos [ cos j+1 cos ( j+1 – 0) / cos max ] . (157) 

 i” = arc cos ( sin i / sin max ) . (158) 

 i” = arc cos [ cos i cos ( i– 0) / cos max ] . (159) 

The next relationship needed for the transformation takes advantage of the fact that the primed 

and double-primed coordinate systems share the same y axis. So in general, y’ = y”. Using the 

definition of equation 2, we may establish that 

 y’ = ( Rave + E ) cos ' sin ' = y” = ( Rave + E ) cos ” sin ”. (160) 

For the case where the point is on the double-primed equator, we may take ” as zero. After 

dividing by ( Rave + E ), we see for a point on the double-primed equator 

 cos ' sin ' = sin ”. (161) 

Again, substituting the values for the jth and j+1th waypoints and the interpolated point from 

equations 142 to 147, and solving for the longitude in the double-primed coordinate system, we 

find that 

 ”j = arc sin [cos j sin ( j – 0 ) ] . (162) 

 ”j+1 = arc sin [cos j+1 sin (  j+1 – 0 ) ] . (163) 

 ”i = arc sin [cos i sin ( i – 0 ) ] . (164) 

Note that equations 162–164 use the inverse sine, or arc sin, function to obtain the longitude in 

the double-primed coordinate system ”, while equations 154–159 use the inverse cosine, or arc 

cos, function to get ”. At first glance, these two functions appear redundant. But the inverse 
68

sine and cosine functions deliver not one but two answers in the range –180° to +180°. This is 

illustrated in figure 45 (28). As an example, dots are placed on the sine and cosine functions 

equal to 0.5, which correspond to an angles of 30° and 150° for sine, and 60° and –60° for cosine. 

Figure 45. The dots on the sine and cosine graph at 0.5 show that there are two distinct angles for any one 

value of sine and cosine in the interval –180° to +180°. 
69

To sort out which of the two angles that each of the inverse functions produce is the correct 

angle, we must first examine the sign of the functions in each of the four quadrants. This is 

shown in figure 46. In the case where the angle  is between 0° and 90° (quadrant I), then cos 

and sin  are both positive. When  is between 90° and 180° (quadrant II), then sin  is positive 

again but cos  is negative. If  is between 0° and –90° (quadrant IV), then cos  is positive but 

sin  is negative. Finally, if  is between –90° and –180° (quadrant III), then both cos  and sin 

 are negative. 

Figure 46. For an angle  ranging from –180° to +180°, the sine and cosine functions take 

on either a positive or negative sign in each 90° quadrant. 

Next, we must consider which values of  are returned by programming languages for the 

functions arc cos and arc sin. As is the case for Java (29), C++ (30, 31), and FORTRAN (the 

Formula Translator) (32),  = arc cos x returns a value for  that is in the first and second 

quadrants, that is 0 > , while  = arc sin x returns a value for that is in the first and 

fourth quadrants, that is –90° <  < +90°. Yet  could be in any of the four quadrants of figure 

46. We see that for all angles of  for which sin  is positive (those in quadrants I and II), the arc 

cos function returns the value of  correctly. For the other two quadrants (III and IV), we see 

that the sine function is negative, but the arc cos function returns a value of that is opposite in 

sign to the values nominally returned using programming languages. Therefore, to make sure the 
70

angle  returned is the correct case, we leave the sign for  as is for the function arc cos x if the 

sin function is positive and change the sign for arc cos X if the sin function is negative. 

Applying equations 154–159 and 162–164, we find that 

 if cos j sin ( j – 0 ) ≥ 0 , 

j” = arc cos ( sin j / sin max ) , 

and 

 if cos j sin ( j – 0 ) < 0 , 

 j” = –arc cos ( sin j / sin max ) . (165) 

Likewise,

 if cos j+1 sin (  j+1 – 0 ) ≥ 0 

 j+1” = arc cos ( sin j+1 / sin max ) , 

and 

 if cos j+1 sin (  j+1 – 0 ) < 0 

 j+1” = –arc cos ( sin j+1 / sin max ) . (166) 

And finally, 

 if cos i sin ( i – 0 ) ≥ 0

 i” = arc cos ( sin i / sin max ) , 

and 

 if cos i sin ( i – 0 ) < 0 

 i” = –arc cos ( sin i / sin max ) . (167) 

We now apply equations 138 and 139 to find the interpolated angle i”. They are repeated here 

as equations 168 and 169. 

 i = Di / (Rave + Ej) . (168) 

 i = [(1 – k2)

1/2 / k] Ln [(Rave + Ei) / (Rave + Ej )] . (169) 

Equation 168 is used if the elevations at the waypoints Ej and Ej+1 are equal. If they are not, then 

equation 169 is used, where the value for k is defined by equation 136, repeated here as equation 

170: 

 k = (Ej+1 – Ej) / Dj . (170) 
71

Notice that it is Dj in equation 170, the distance between the waypoints, not to be confused with 

the interpolated distance Di, the distance from the jth waypoint to the platform at the interpolated 

time ti. One more item needed is the angular distance between the waypoints. That is, the angle 

between two lines: one that goes from waypoint j to the Earth’s center, and the other that goes 

between the j+1th waypoint to the Earth’s center. This angle, , was calculated in equation 80 

and is repeated here as equation 171. 

  = 2 arc sin { [( 1 / 2 ) [(cos j+1 cos j+1 – cosj cos j)

2

 + 

 (cos j+1 sin j+1 – cos j sin j)

2

 + (sin j+1 – sin j)

2

]

1/2} . (171) 

(j,j) and (j+1,j+1) are the longitude and latitude of the jth and j+1th waypoints in the original 

coordinate system. The interpolated angle in the double-primed coordinate system ” is the 

fraction of the angular distance  between the waypoints that the platform reached at the 

interpolated time ti. Expressed quantitatively using the terms derived so far, it is the fraction 

i / . The radial distance covered between the waypoints by the platform is the fraction of 

the distance covered times the difference between the coordinates in the double-primed 

coordinate system, or ( j+1”– j” ) ( i /  ). To obtain the double-primed interpolated 

longitude, we need only add this expression to the longitude of the jth waypoint. 

 i” = j” + ( j+1”– j” ) ( i /  ) . (172) 

Note that the quantities i and  are always positive. But depending on the direction of travel, 

j” may be greater than or less than j+1”. To make sure that i” is between j” and j+1”, it is 

necessary to use the full expression in equation 172 rather than simply add i to j”. As a final 

step, we make sure that the interpolated longitude in the double-primed coordinate system i” is 

between –180° and +180°. If it is not, then 360° is either added or subtracted to the value of i” 

until it is. 

Now that we have the interpolated point (,i”) in the double-primed coordinate system 

(remember, the point is on the double-primed equator), we reverse the two-step process to obtain 

the coordinates (i,i) in the original coordinate system. We shall review equations 158, 164, 

and 155. Solving equation 158 for i, we see that 

 i = arc sin (cos i” sin max) . (173) 

The range of i is between –90° and +90°, and we see that equation 173 delivers a unique value 

for i. Now that we know the value of i, we may solve for i. Solving equations 164 and 155 

for I, we find 

 i = arc sin (sin i” / cos i) + 0 (174) 

and 

 i = arc cos (cos i” cos max / cos i ) + 0 . (175)
72

Once again, we find that i ranges from –180° to + 180°. But the range of values offered by 

computer languages (29–32) returns values between –90° and +90° for arc sin and 0° and 180° 

for arc cos. Reusing the process when we transformed to the double-primed coordinate system, 

we arrive at the value of i: 

 if sin i” / cos i ≥ 0, then 

 i = arc cos (cos i” cos max / cos i ) + 0 , 

 and if sin i” / cos i < 0, then 

 i = –arc cos (cos i” cos max / cos i ) + 0 . (176) 

Keep in mind that i must be between –180° and +180°. If it is not, then add or subtract 360° as 

is appropriate. We have finished describing the process of obtaining interpolated points along a 

great circle on a sphere. 

4.3.5 Summary of Inputting Waypoint Data and Obtaining Interpolated Points 

This section briefly describes the process of inputting waypoint data and obtaining the location 

of the interpolated points along a great circle. 

Step 1: Assign the waypoint longitudes j, latitudes j, and elevations Ej. If the waypoints are 

for a ground platform, the elevations are read from a database. 

Step 2: Calculate the great circle equation that connects consecutive waypoints. Most of the 

time, it will be possible to use equation 10 to get 0 and equation 11 to get max, which are the 

descriptive parameters of equation 6, the great circle equation. 

Step 3: Calculate the direct distance with equation 79, the angular distance  with equation 

80, and the constant k with equation 95 for consecutive waypoints. This makes it possible to 

calculate the distance Dj between consecutive waypoints with equation 98. 

Step 4: Obtain the cruising velocity vcru and acceleration constant a0 for the platform from a 

database. Then, assign either the arrival time ta,j+1 to each waypoint to calculate the waypoint 

velocity vj+1 (which can usually be done with equation 76) or assign the waypoint velocity and 

calculate the arrival time ta,j+1 (which can usually be done with equation 75). 

Step 5: Calculate the inflection times t1 and t2, which can usually be found using equations 72 

and 73. 

Step 6: Select the loiter time (td,j – ta,j) for each waypoint. Notice that this can only be done if 

the waypoint velocity vj is zero. If the waypoint velocity vj is not zero, the loiter time defaults to 

zero.
73

Note that steps 1–6 are executed every time a waypoint is added to a platform path. After all the 

waypoints for a platform have been entered, proceed with step 7–10 and calculate the 

interpolated latitudes and longitudes (i,i). 

Step 7: Select the interpolation times, (each one delineated ti) for the platform. Use these times 

to calculate the interpolated distances Di with equations 133–135. Obtain the interpolated 

elevation Ei with equation 137 and the interpolated angle i with equation 139. 

Step 8: Transform the waypoints with departure times td,j and arrival times ta,j+1, which bound 

the interpolated time ti to a new, double-primed coordinate system with longitudes j” and j+1” 

using equations 165 and 166. The points will be on the double-primed equator. 

Step 9: Calculate the interpolated longitude i” using equation 172. 

Step 10: Transform the interpolated longitude i” (it too is on the double-primed equator) back 

to the original coordinate system with the longitude and latitude (i,i) using equations 173 and 

176. 

5. Conclusions 

The algorithms and formulae describing the motion of communications platforms needed for the 

NCAM deployment module have been successfully derived. For purposes of verification, the 

details of the algorithm and formulae development have been documented and are shown to be 

consistent with the assumptions implicit in the model’s description of platform motion. The 

assumptions were that the platform travels along great circle paths between waypoints, the 

platforms accelerate and decelerate at a constant rate between waypoints, the platforms have a 

maximum cruising speed consistent with the platform and terrain type (HMMWV vs. tank vs. 

UAV vs. dismounted Soldier on paved road vs. rough mountainous terrain vs. flat desert terrain, 

etc.), and the platforms change elevation at a constant rate. 

This approach allowed the user to outline each platform’s path with waypoints consistent with 

the user’s specification of arrival and departure times at each waypoint. Using the equations 

derived provided the user with upper and lower limits for inputted data. This gave the user 

maximum freedom to select input and ensure that the inputs were restricted so as to be selfconsistent. This made it possible to calculate the position of platforms between waypoints at 

predetermined snapshot times. Hence, it will become possible to model link viability between 

moving platforms once the other modules are completed. Therefore, the equations and 

algorithms presented here are mandatory to allowing NCAM to model communications links 

between moving nodes. 
74

Making the motion model more complex would have made it possible to model platform motion 

more realistically. However, the derivations of the equations and algorithms would have to be 

far more involved and complex than the ones shown here. Furthermore, it would have demanded 

that more information be input from the user. Making the motion simpler would have 

diminished the realism of the platform motion and diminished the control the user has over 

platform motion. The approach used here (specifying waypoint location and the 

arrival/departure velocity or time) provided the best compromise between user control, equation 

and algorithm simplicity, and realistic platform motion.
75

6. References 

1. Still, G. W.; Nealon, J. F. The Case for Using the Spherical Model to Calculate the 

Interpolated Points in the Connectivity Software Deployment Module; ARL-TR-4373; U.S. 

Army Research Laboratory: Aberdeen Proving Ground, MD, February 2008. 

2. Eppink, D.; Keubler, W.; TIREM/SEM Handbook; ECAC-HDBK-93-076; Department of 

Defense Electromagnetic Compatibility Analysis Center: Annapolis, MD, March 1994; pp 

1–134. 

3. Wikipedia. http://en.wikipedia.org/wiki/Great_circle (accessed 25 October 2007), search: 

Great Circle. 

4. Math Is Fun. http://www.mathisfun.com/equation_of_line.html (accessed 21 September 

2007), search: Equation of a Straight Line. 

5. S.O.S. Mathematics. http://www.sosmath.com/trig/Trig5/trig5/trig5.html (accessed 17 

October 2007), search: Trig Homework Solver. 

6. Wikipedia. http://en.wikipedia.org/wiki/Trigonometric_function (accessed 4 April 2008), 

search: Trigonometric Functions. 

7. Cress, P.; Dirksen, P.; Graham, J. W. Fortran IV With Watfor and Watfiv; Prentice-Hall 

Inc.: Englewood Cliffs, NJ, 1970; p 407. 

8. Hofstra People. http://people.hofstra.edu/Stefan_Waner/trig/trig2.html (accessed 9 April 

2008), search: The Trigonometric Functions. 

9. RCN Corporation Home Page. http://users.rcn.com/mwhitney.massed/trigresources/ 

trig-graph-properties.html (accessed 16 January 2008), search: Properties of the Basic 

Trigonometric Functions. 

10. Analyze Math. http://www.analyzemath.com/trigonometry/properties.html (accessed 16 

January 2008), search: Properties of Trigonometric Functions. 

11. Wikipedia. http://en.wikipedia.org/wiki/Oblate_spheroid (accessed 7 February 2008), 

search: Oblate Spheroid. 

12. Wikipedia. http://en.wikipedia.org/wiki/Ellipsoid (accessed 7 February 2008), search: 

Ellipsoid. 

13. Wikipedia. http://en.wikipedia.org/wiki/WGS84 (accessed 7 February 2008), search: World 

Geodetic System.
76

14. Ed Williams’ Aviation Home Page. http://williams.best.vwh.net/ellipsoid/node1.html 

(accessed 10 December 2007), search: Spheroid Geometry. 

15. Analyze Math. http://www.analyzemath.com/EllipseEq/EllipseEq.html (accessed 10 

December 2007), search: Equation of Ellipse. 

16. National Geodetic Survey Home Page. http://www.ngs.noaa.gov/PC_PROD/Inv_Fwd/ 

(accessed 29 June 2007), search: PC Software Download-INVERSE and FORWARD. 

17. Everything 2. http://everything2.com/index.pl?node_id=1056913 (accessed 20 February 

2008), search: Differential Solid Angle. 

18. Math World. http://mathworld.wolfram.com/OblateSpheroid.html (accessed 20 February 

2008), search: Oblate Spheroid. 

19. Grobner, W.; Hofreiter, N. Integraltafel erster teil Unbestimmte Integrale, No. 5a; 

Springer-Verlang: New York, 1975; p 45. 

20. Dwight, H. W. Tables of Integrals and Other Mathematical Data, 4th ed., No. 601.2; 

MacMillan Publishing Company, Inc.: New York, 1961; p 137. 

21. Weast, R., Astile, M., Eds.; CRC Handbook of Chemistry and Physics, 63rd ed.; CRC Press 

Inc.: Boca Raton, FL, 1983; p F-155. 

22. Mission-Planning. http://www.mission-planning.com/DTED_Part2.htm (accessed 28 June 

2007), search: Pablo’s Mission Planning Website. 

23. Wolfram Demonstrations Project Home Page. http://demonstrations.wolfram.com/ 

ShortestPathBetweenTwoPointsOnASphere/ (accessed 27 February 2008), search: Shortest 

Path Between Two Points on a Sphere. 

24. University of California Santa Cruz Mathematics Home Page. http://math.ucsc.edu 

/~bcauchon/math22_eqns.pdf (accessed 21 March 2008), search: Math 22 Equations. 

25. Leithold, L. The Calculus With Analytic Geometry, 3rd ed.; Harper & Row Publishers: New 

York, NY, 1976; p 784. 

26. Wikipedia. http://en.wikipedia.org/wiki/Natural_logarithm (accessed 9 April 2008), search: 

Natural Logarithm. 

27. The Math Page. http://www.themathpage.com/aTrig/arc-length.htm (accessed 10 April 

2008), search: The Math Page Trigonometry 

28. Wikipedia. http://en.wikipedia.org/wiki/Image:Sine_cosine_plot.svg (accessed 18 

September 2008), search: Sine Cosine Plot.sfg. 

29. Sun Spot World. http://www.sunspotworld.com/docs/Purple/javadoc/com/sun/ 

squawk/util/MathUtils.html (accessed 18 September 2008), search: Class Math Utils. 
77

30. Salt Lake City School District Home Page. http://www.slc.k12.ut.us/staff/perickson 

/programming/tutorials/acos.html (accessed 18 September 2008), search: acos. 

31. Salt Lake City School District Home Page. http://www.slc.k12.ut.us/staff/perickson 

/programming/tutorials/asin.html (accessed 18 September 2008), search: asin. 

32. Cress, P.; Dirksen, P.; Graham, J. W. Fortran IV With Watfor and Watfiv; Prentice Hall, 

Inc.: Englewood Cliffs, NJ, 1970; p 407. 
78

INTENTIONALLY LEFT BLANK. 
79

List of Symbols, Abbreviations, and Acronyms 

CCI connectivity confidence interval 

DTED Digital Terrain Elevation Data 

FORTRAN Formula Translator 

HMMWV high-mobility multipurpose wheeled vehicle 

NCAM Network Connectivity Analysis Model 

S/N signal to noise 

TIREM Terrain Integrated Rough Earth Model 

UAV unmanned aerial vehicle 

WGS84 World Geodetic System for 1984
NO. OF 

COPIES ORGANIZATION 

80

 1 DEFENSE TECHNICAL 

 (PDF INFORMATION CTR 

 only) DTIC OCA 

 8725 JOHN J KINGMAN RD 

 STE 0944 

 FORT BELVOIR VA 22060-6218 

 1 DIRECTOR 

 US ARMY RESEARCH LAB 

 IMNE ALC HRR 

 2800 POWDER MILL RD 

 ADELPHI MD 20783-1197 

 1 DIRECTOR 

 US ARMY RESEARCH LAB 

 RDRL CIM L 

 2800 POWDER MILL RD 

 ADELPHI MD 20783-1197 

 1 DIRECTOR 

 US ARMY RESEARCH LAB 

 RDRL CIM P 

 2800 POWDER MILL RD 

 ADELPHI MD 20783-1197 

 1 DIRECTOR 

 US ARMY RESEARCH LAB 

 RDRL D 

 2800 POWDER MILL RD 

 ADELPHI MD 20783-1197 

ABERDEEN PROVING GROUND

 1 DIR USARL 

 RDRL CIM G (BLDG 4600) 
NO. OF 

COPIES ORGANIZATION 

81

 1 ASST SECY ARMY 

 (CD ACQSTN LOGISTICS & TECH 

 only) SAAL ZP RM 2E661 

 103 ARMY PENTAGON 

 WASHINGTON DC 20310-0103 

 1 ASST SECY ARMY 

 (CD ACQSTN LOGISTICS & TECH 

 only) SAAL ZS RM 3E448 

 103 ARMY PENTAGON 

 WASHINGTON DC 20310-0103 

 1 DIRECTOR FORCE DEV 

 (CD DAPR FDZ 

 only) RM 3A522 

 460 ARMY PENTAGON 

 WASHINGTON DC 20310-0460 

 1 US ARMY TRADOC ANL CTR 

 ATRC W 

 A KEINTZ 

 WSMR NM 88002-5502 

 1 USARL 

 RDRL SLE 

 R FLORES 

 WSMR NM 88002-5513 

ABERDEEN PROVING GROUND

 1 US ARMY DEV TEST COM 

 CSTE DTC TT T 

 314 LONGS CORNER RD 

 APG MD 21005-5055 

 1 US ARMY EVALUATION CTR 

 CSTE AEC SVE 

 R LAUGHMAN 

 4120 SUSQUEHANNA AVE 

 APG MD 21005-3013 

 12 DIR USARL 

 RDRL SL 

 J BEILFUSS 

 J FEENEY 

 J FRANZ 

 M STARKS 

 P TANENBAUM 

 RDRL SLB A 

 G MANNIX 

 RDRL SLB A 

 D FARENWALD 

 RDRL SLB D 

 R GROTE 

 RDRL SLB E 

 M PERRY 

 RDRL SLB G 

 P MERGLER 

 RDRL SLB S 

 S SNEAD 

 RDRL SLB W 

 L ROACH 
NO. OF 

COPIES ORGANIZATION 

82

 1 US ARMY RSRCH DEV & ENGRG 

CMD 

 SYS OF INTEGRATION 

 AMSRD SS T 

 6000 6TH ST STE 100 

 FORT BELVOIR VA 22060-5608 

 1 US ARMY RSCH LAB 

 RDRL SLE G 

 L ANDERSON 

 WSMR NM 88002-5502 

 1 US ARMY RSCH LAB 

 RDRL SLE M 

 P SIMPSON 

 WSMR NM 88002-5502 

 5 US ARMY RSCH LAB 

 RDRL SLE W 

 P BOTHNER 

 FORT MONMOUTH NJ 07703-5602 

 1 US ARMY RDECOM TARDEC 

 S CAMPBELL 

 MS 263 

 6501 E 11 MILE RD 

 WARREN MI 48397-5000 

ABERDEEN PROVING GROUND

 6 DIR USARL 

 RDRL SL 

 D BAYLOR 

 RDRL SLE M 

 G W STILL (5 CPS)