geodesics  R Documentation 
Function geodesics
carries out the calculation of pairwise geodesic
distances within a set of coordinates or between two sets thereof, using one
of two calculation approaches.
geodesics(
x,
y,
method = c("haversine", "Vincenty"),
radius = 6371000,
sma = 6378137,
flat = 1/298.257223563,
maxiter = 1024L,
tol = .Machine$double.eps^0.75
)
x 
A set of geographic coordinates in the form of a twocolumn

y 
An other twocolumn 
method 
The calculation method used to obtain the distances (default: haversine method; see details). 
radius 
Radius of the planetary body (when assuming a sphere; default: 6371000 m). 
sma 
Length of the semimajor axis of the planetary body (when assuming a revolution ellipsoid; default: 6378137 m). 
flat 
Flattening of the ellipsoid (default: 1/298.257223563). 
maxiter 
Maximum number of iterations, whenever iterative calculation is involved (default: 1024). 
tol 
Tolerance used when iterative calculation is involved (default:

When only one set of coordinates is given to the function (i.e.,
when argument y
is omitted), the function returns the pairwise
distances in the form of a 'distclass' object representing a
lowertriangle matrix. When the second coordinate set is given, the function
calculates the distances between each coordinate of argument x
and
each coordinate of argument y
.
Two calculation methods are implemented. The first is the haversine formula,
which assume the planetary body to be a sphere. The radius of that sphere is
given to the function as its argument radius
, with the default value
being the mean radius of planet earth. Of the two methods implemented, the
haversine formula is fastest but its precision depend on how well the
planetary body match the sphericity assumption. The second method implemented
is Vincenty's inverse formula, which assumes the the planetary body is a
revolution ellipsoid, which is expected for rotating semifluid such as
planet earth. Argument sma
, the length of the semimajor axis,
corresponds to the radius of the circle obtained when the revolution
ellipsoid at the equator, whereas argument flat
correspond to the
compression of the sphere, along the diameter joining the poles, to form the
ellipsoid of revolution. Their default values corresponds to parameters for
planet Earth according to WGS84. These values, along with arguments
maxiter
and tol
, are ignored when using the haversine formula,
whereas the value of argument radius
is ignored when using Vincenty's
inverse formula.
Vincenty's inverse formula is more precise on planet Earth (on the order of
0.5mm) than the haversine formula, but it involves more computation time and
may sometimes fail to converge. This is more likely for pairs of locations
that are nearly antipodal or both (numerically) very close to the equator.
The results returned by the function when using Vincenty's inverse formula
are given a niter
attribute that gives the number of iterations that
were necessary to achieve convergence. Numbers greater than argument
maxiter
are indicative of failed convergence; a warning is issued in
such a circumstance.
Geodesic distance matrices are non metric.
A 'distclass' object or, whenever argument y
is
provided, a matrix
with as many rows as the number of rows in
argument x
and as many columns as the number of rows in argument
y
.
Guillaume Guenard and Pierre Legendre, Bertrand Pages Maintainer: Guillaume Guenard <guillaume.guenard@gmail.com>
Vincenty, T. 1975. Direct and Inverse Solutions of Geodesics on the Ellipsoid with application of nested equations. Survey Review XXIII (176): 8893 doi:10.1179/sre.1975.23.176.88
Inman, J. 1835. Navigation and Nautical Astronomy: For the Use of British Seamen (3 ed.). London, UK: W. Woodward, C. & J. Rivington
The distclass
and associated methods.
##
### First example: locations spread throughout the world
##
coords < cbind(c(43,22,9,12,40,72,86,22),
c(135,22,0,1,45,12,27,139))
res_hav < geodesics(coords) ## Default: the haversine formula
res_hav
res_vif < geodesics(coords, method = "Vincenty")
res_vif
attr(res_vif,"niter") ## The numbers of iterations
res_vifres_hav ## Absolute difference
200*(res_vifres_hav)/(res_vif+res_hav) ## Large relative difference
##
### Second example: locations nearer from one another
##
coords < cbind(c(45.01,44.82,45.23,44.74),
c(72.03,72.34,71.89,72.45))
res_hav < geodesics(coords)
res_vif < geodesics(coords, method = "Vincenty")
res_vifres_hav ## Absolute difference
200*(res_vifres_hav)/(res_vif+res_hav) ## Relative difference are smaller
##
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.