如何计算R中两个纬度/经度之间的方位

Question

我想计算两个位置之间的方位(方向(。这就是样本的样子：

# A tibble: 10 x 14
   vehicle_id                  trip_id seg_1_lat seg_1_lon seg_2_lat seg_2_lon seg_3_lat seg_3_lon seg_4_lat seg_4_lon seg_5_lat seg_5_lon seg_6_lat seg_6_lon
   <chr>                         <int>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>
 1 Zz3yE90z++QmTX2QO5dHI78IK/Q       1      13.7     101.       13.7     101.       13.6     101.       13.7     101.       13.7     101.       13.7     101. 
 2 Zz3yE90z++QmTX2QO5dHI78IK/Q       2      13.7     101.       13.7     101.       13.7     101.       13.7     101.       13.7     101.       13.7     101. 
 3 Zz3yE90z++QmTX2QO5dHI78IK/Q       3      13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101. 
 4 Zz3yE90z++QmTX2QO5dHI78IK/Q       4      13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101. 
 5 Zz3yE90z++QmTX2QO5dHI78IK/Q       5      13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101. 
 6 Zz3yE90z++QmTX2QO5dHI78IK/Q       6      13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101. 
 7 Zz3yE90z++QmTX2QO5dHI78IK/Q       7      13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101.       13.6     101. 
 8 /+bx80f3gOoPMoFBsS+3xX6jpi8   44496     101.       13.9     101.       13.9     101.       13.9     101.       13.9     101.       14.0     101.       14.0
 9 /+bx80f3gOoPMoFBsS+3xX6jpi8   44497     101.       14.0     101.       14.0     101.       14.0     101.       14.0     101.       14.0     101.       14.0
10 /+bx80f3gOoPMoFBsS+3xX6jpi8   44498     101.       13.8     101.       13.8     101.       13.9     101.       13.9     101.       13.9     101.       13.9

my.df <- structure(list(vehicle_id = c("Zz3yE90z++QmTX2QO5dHI78IK/Q","Zz3yE90z++QmTX2QO5dHI78IK/Q", "Zz3yE90z++QmTX2QO5dHI78IK/Q","Zz3yE90z++QmTX2QO5dHI78IK/Q", "Zz3yE90z++QmTX2QO5dHI78IK/Q","Zz3yE90z++QmTX2QO5dHI78IK/Q", "Zz3yE90z++QmTX2QO5dHI78IK/Q","/+bx80f3gOoPMoFBsS+3xX6jpi8", "/+bx80f3gOoPMoFBsS+3xX6jpi8","/+bx80f3gOoPMoFBsS+3xX6jpi8"), trip_id = c(1L, 2L, 3L, 4L, 5L,6L, 7L, 44496L, 44497L, 44498L), seg_1_lat_lat = c(13.67829,13.66859, 13.63345, 13.60629, 13.56958, 13.60148, 13.60283, 13.93405,14.00297, 13.80508), seg_1_lon_lon = c(100.62387, 100.72032,100.71102, 100.69999, 100.78749, 100.71103, 100.70543, 100.72371,100.6638, 100.71231), seg_2_lat_lat = c(13.65382, 13.66353, 13.63349,13.61857, 13.57759, 13.60519, 13.60118, 13.92169, 13.99238, 13.82186), seg_2_lon_lon = c(100.67562, 100.71894, 100.71096, 100.70363,100.78828, 100.7129, 100.71196, 100.69827, 100.65185, 100.72064), seg_3_lat_lat = c(13.63679, 13.66353, 13.63349, 13.6344, 13.58393,13.61076, 13.60193, 13.92155, 13.98482, 13.85196), seg_3_lon_lon = c(100.71057,100.71894, 100.71096, 100.71093, 100.79077, 100.71317, 100.71915,100.68866, 100.61921, 100.72147), seg_4_lat_lat = c(13.65828,13.66353, 13.63349, 13.6344, 13.58385, 13.62504, 13.60487, 13.92121,13.98708, 13.88833), seg_4_lon_lon = c(100.71631, 100.71894,100.71096, 100.71093, 100.79403, 100.71832, 100.71483, 100.65924,100.61634, 100.72214), seg_5_lat_lat = c(13.65828, 13.65958,13.63349, 13.6344, 13.59328, 13.63082, 13.605, 13.95379, 13.98707,13.92194), seg_5_lon_lon = c(100.71631, 100.71773, 100.71096,100.71093, 100.76465, 100.71459, 100.71402, 100.62163, 100.61612,100.70304), seg_6_lat_lat = c(13.65828, 13.65258, 13.63349, 13.63433,13.59612, 13.63241, 13.60567, 13.95905, 13.98712, 13.90366),seg_6_lon_lon = c(100.71631, 100.71571, 100.71096, 100.71101,100.74922, 100.71297, 100.71292, 100.62037, 100.61631, 100.66788)), row.names = c(NA, -10L), class = c("tbl_df", "tbl", "data.frame"))

因此，它将在两个位置之间进行计算。例如，seg_1 with seg_2和seg_2 with seg_3以及seg_3 with seg_4和seg_4 with seg_5以及seg_5 with seg_6。因此，将产生5个方向的结果。

我期望的输出应该有这样的列；

vehicle_id                     trip_id    bearing_1   bearing_2  bearing_3  bearing_4  bearing_5                                        
<chr>                         <dbl>       <dbl>       <dbl>      <dbl>      <dbl>      <dbl>

Answer 1

使用for循环和包geosphere可以执行：

for (i in 1:5) 
     my.df[[paste0("bearing_", i)]] <- 
         bearing(my.df[,paste0("seg_",i,c("_lon_lon","_lat_lat"))], 
                 my.df[,paste0("seg_",i+1,c("_lon_lon","_lat_lat"))])

将列添加到my.df。请注意，有些行的连续值是相同的-方位未定义，但geosphere返回-180。如果这对你的分析很重要，你应该明确地检查这些。

Answer 2

试试这个：

out <- as.data.frame(lapply(c(3,5,7,9,11),
   function(cn) geosphere::bearing(as.matrix(my.df[,cn+c(1,0)]), 
                                   as.matrix(my.df[,cn+c(3,2)]))))
names(out) <- paste0("bearing_", seq_along(out))
cbind(my.df[,1:2], out)
#                     vehicle_id trip_id   bearing_1   bearing_2   bearing_3  bearing_4  bearing_5
# 1  Zz3yE90z++QmTX2QO5dHI78IK/Q       1  115.799976  116.480619   14.638536 -180.00000 -180.00000
# 2  Zz3yE90z++QmTX2QO5dHI78IK/Q       2 -165.067040 -180.000000 -180.000000 -163.32392 -164.24079
# 3  Zz3yE90z++QmTX2QO5dHI78IK/Q       3  -55.719426 -180.000000 -180.000000 -180.00000 -180.00000
# 4  Zz3yE90z++QmTX2QO5dHI78IK/Q       4   16.167746   24.275570 -180.000000 -180.00000  131.81829
# 5  Zz3yE90z++QmTX2QO5dHI78IK/Q       5    5.510892   21.016200   91.436657  -71.83082  -79.34176
# 6  Zz3yE90z++QmTX2QO5dHI78IK/Q       6   26.244111    2.714509   19.429022  -32.25773  -44.89851
# 7  Zz3yE90z++QmTX2QO5dHI78IK/Q       7  104.483610   83.911999  -55.170513  -80.68170  -58.08888
# 8  /+bx80f3gOoPMoFBsS+3xX6jpi8   44496 -116.443052  -90.853309  -90.674317  -48.42547  -13.16726
# 9  /+bx80f3gOoPMoFBsS+3xX6jpi8   44497 -132.224000 -103.339530  -51.117070  -92.66506   74.91771
# 10 /+bx80f3gOoPMoFBsS+3xX6jpi8   44498   25.878444    1.543343    1.031031  -29.03379 -118.01960