Great circle distance matrix or vector
Given two sets of longitude/latitude locations, rdist.earth
computes
the Great circle (geographic) distance matrix among all pairings and
rdist.earth.vec
computes a vector of pairwise great circle distances
between corresponding elements of the input locations using the Haversine
method and is used in empirical variogram calculations.
rdist.earth(x1, x2, miles = TRUE, R = NULL) RdistEarth(x1, x2=NULL, miles=TRUE, R=NULL) rdist.earth.vec(x1, x2, miles = TRUE, R = NULL)
x1 |
Matrix of first set of lon/lat coordinates first column is the longitudes and second is the latitudes. |
x2 |
Matrix of second set of lon/lat coordinates first column is the longitudes and second is the latitudes. If missing or NULL x1 is used. |
miles |
If true distances are in statute miles if false distances in kilometers. |
R |
Radius to use for sphere to find spherical distances. If NULL the radius is either in miles or kilometers depending on the values of the miles argument. If R=1 then distances are of course in radians. |
Surprisingly the distance matrix is computed efficiently in R by dot products of the
direction cosines. This is the calculation in rdist.earth
. Thanks to Qing Yang for pointing this out a long time
ago. A more efficient version has been implemented in C with the
R function RdistEarth
by Florian Gerber who has also experimented with parallel versions of fields functions.
The main advantage of RdistEarth
is the largely reduce memory usage.
The speed seems simillar to rdist.earth
. As Florian writes:
"The current fields::rdist.earth() is surprisingly fast. In the case where only the argument 'x1' is specified, the new C implementation is faster. In the case where 'x1' and 'x2' are given, fields::rdist.earth() is a bit faster. This might be because fields::rdist.earth() does not check its input arguments and uses a less complicated (probably numerically less stable) formula."
The great circle distance matrix if nrow(x1)=m and nrow( x2)=n then the returned matrix will be mXn.
Doug Nychka, John Paige, Florian Gerber
data(ozone2) out<- rdist.earth ( ozone2$lon.lat) #out is a 153X153 distance matrix out2<- RdistEarth ( ozone2$lon.lat) all.equal(out, out2) upper<- col(out)> row( out) # histogram of all pairwise distances. hist( out[upper]) #get pairwise distances between first 10 and second 10 lon/lat points x1 = ozone2$lon.lat[1:10,] x2 = ozone2$lon.lat[11:20,] dists = rdist.earth.vec(x1, x2) print(dists)
Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.