1

我确实有带有 gps 数据的数据框，我需要计算两行之间的距离（当前和上一个）

            id                 time       lat      long heartrate altitude
1  20130424.tcx 2013-04-24T04:53:22Z 50.024818 14.522724       151     <NA>
2  20130424.tcx 2013-04-24T04:53:26Z 50.024818 14.522724        96     <NA>
3  20130424.tcx 2013-04-24T04:53:30Z 50.024818 14.522724       104     <NA>
4  20130424.tcx 2013-04-24T04:53:34Z 50.024818 14.522724       107     <NA>
5  20130424.tcx 2013-04-24T04:53:38Z 50.024818 14.522724       108     <NA>
6  20130424.tcx 2013-04-24T04:53:42Z 50.024818 14.522724       112    372.0
7  20130424.tcx 2013-04-24T04:53:46Z 50.024818 14.522724       151    372.0
8  20130424.tcx 2013-04-24T04:53:47Z 50.024677 14.522874       151    356.0
9  20130424.tcx 2013-04-24T04:53:50Z 50.024677 14.522874       118    356.0
10 20130424.tcx 2013-04-24T04:53:54Z 50.024677 14.522874       118    356.0
11 20130424.tcx 2013-04-24T04:53:58Z 50.024464 14.522917       147    358.0
12 20130424.tcx 2013-04-24T04:54:02Z 50.024464 14.522917       144    358.0
13 20130424.tcx 2013-04-24T04:54:06Z 50.024269 14.522853       150    367.0
14 20130424.tcx 2013-04-24T04:54:10Z 50.024269 14.522853       152    367.0
15 20130424.tcx 2013-04-24T04:54:13Z 50.024002 14.522874       152    380.0

我能够将数据连接到自身并获取每一行的前一行（可能有更简单的解决方案）：

library(sqldf)
mydft    = mydf[-nrow(mydf),] 
mydft$id  = seq_along(mydft$id) +1
mydf$id  = seq_along(mydf$id)
mydft2 <- sqldf("select a.*, b.lat as lat2, b.long as long2 from mydf a left join mydft b using (id)")

我现在如何计算与列的距离lat, long, lat2, long2？我尝试过这里描述的方法：

R <- 6371 # Earth mean radius [km]
mydft2$delta.long <- (mydft2$long2 - mydft2$long)
mydft2$delta.lat <- (mydft2$lat2 - mydft2$lat)
mydft2$a <- sin(mydft2$delta.lat/2)^2 + cos(mydft2$lat) * cos(mydft2$lat2) * sin(mydft2$delta.long/2)^2
mydft2$c <- 2 * asin(min(1,sqrt(mydft2$a)))
mydft2$d = R * c

但这仅返回列表 NA。

4

1 回答 1

5

一种简单地实现这一目标的方法是使用R设施进行空间分析。对于这个例子，我们可以使用spDistN1优秀sp包中的函数。

第一步是将您的数据转换为 SpatialPoints（或 SpatialPointsDataFrame），我在此示例中假设您的点位于地理投影中（带 WSG84 基准的 longlat）

txt <- "            id                 time       lat      long heartrate altitude
20130424.tcx 2013-04-24T04:53:22Z 50.024818 14.522724       151     <NA>
20130424.tcx 2013-04-24T04:53:26Z 50.024818 14.522724        96     <NA>
20130424.tcx 2013-04-24T04:53:30Z 50.024818 14.522724       104     <NA>
20130424.tcx 2013-04-24T04:53:34Z 50.024818 14.522724       107     <NA>
20130424.tcx 2013-04-24T04:53:38Z 50.024818 14.522724       108     <NA>
20130424.tcx 2013-04-24T04:53:42Z 50.024818 14.522724       112    372.0
20130424.tcx 2013-04-24T04:53:46Z 50.024818 14.522724       151    372.0
20130424.tcx 2013-04-24T04:53:47Z 50.024677 14.522874       151    356.0
20130424.tcx 2013-04-24T04:53:50Z 50.024677 14.522874       118    356.0
20130424.tcx 2013-04-24T04:53:54Z 50.024677 14.522874       118    356.0
20130424.tcx 2013-04-24T04:53:58Z 50.024464 14.522917       147    358.0
20130424.tcx 2013-04-24T04:54:02Z 50.024464 14.522917       144    358.0
20130424.tcx 2013-04-24T04:54:06Z 50.024269 14.522853       150    367.0
20130424.tcx 2013-04-24T04:54:10Z 50.024269 14.522853       152    367.0
20130424.tcx 2013-04-24T04:54:13Z 50.024002 14.522874       152    380.0"

gpsdat <- read.table(text = txt, header = TRUE, na.strings = "<NA>")
str(gpsdat)

require(sp)
coordinates(gpsdat) <- ~ long + lat
proj4string(gpsdat) <- CRS("+proj=longlat +datum=WGS84")

在第二步中，我们现在可以使用该spDistN1功能

sapply(seq_along(gpsdat[-1, ]), function(i)
       spDistsN1(pts = gpsdat[i, ], pt = gpsdat[i+1, ], longlat = TRUE))

 [1] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
 [7] 0.019004 0.000000 0.000000 0.023875 0.000000 0.022155
[13] 0.000000 0.029716

根据我们使用的 GPS 数据类型，您可以R使用readOGR函数（包rgdal）直接读取这些类型的数据

于 2013-04-24T15:13:52.313 回答