R version 3.2.5 (2016-04-14)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 7 x64 (build 7601) Service Pack 1
locale:
[1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United
States.1252
[3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C
[5] LC_TIME=English_United States.1252
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] RPyGeo_0.9-3 RSAGA_0.94-5 plyr_1.8.3
gstat_1.1-3
[5] shapefiles_0.7 foreign_0.8-66 rgdal_1.1-8
maptools_0.8-39
[9] geosphere_1.5-1 ggplot2_2.1.0 raster_2.5-2
rgeos_0.3-19
[13] sp_1.2-3 RevoUtilsMath_3.2.5
loaded via a namespace (and not attached):
[1] Rcpp_0.12.4 maps_3.1.0 munsell_0.4.3 colorspace_1.2-6
lattice_0.20-33
[6] FNN_1.1 xts_0.9-7 tools_3.2.5 parallel_3.2.5
grid_3.2.5
[11] gtable_0.2.0 intervals_0.15.1 digest_0.6.9 mapproj_1.2-4
labeling_0.3
[16] scales_0.4.0 spacetime_1.1-5 zoo_1.7-12
Spatial line: l1
str(l1)
Formal class 'SpatialLinesDataFrame' [package "sp"] with 4 slots
..@ data :'data.frame': 1 obs. of 2 variables:
.. ..$ OBJECTID : int 1
.. ..$ SHAPE_Leng: num 1043
.. ..- attr(*, "data_types")= chr [1:2] "N" "F"
..@ lines :List of 1
.. ..$ :Formal class 'Lines' [package "sp"] with 2 slots
.. .. .. ..@ Lines:List of 1
.. .. .. .. ..$ :Formal class 'Line' [package "sp"] with 1 slot
.. .. .. .. .. .. ..@ coords: num [1:2, 1:2] 642236 643006 6197243 6197948
.. .. .. ..@ ID : chr "0"
..@ bbox : num [1:2, 1:2] 642236 6197243 643006 6197948
.. ..- attr(*, "dimnames")=List of 2
.. .. ..$ : chr [1:2] "x" "y"
.. .. ..$ : chr [1:2] "min" "max"
..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
.. .. ..@ projargs: chr "+proj=utm +zone=20 +south +datum=WGS84 +units=m
+no_defs +ellps=WGS84 +towgs84=0,0,0"
Spatial line in geographic coordinates: l1.geo
str(l1.geo)
Formal class 'SpatialLinesDataFrame' [package "sp"] with 4 slots
..@ data :'data.frame': 1 obs. of 2 variables:
.. ..$ OBJECTID : int 1
.. ..$ SHAPE_Leng: num 1043
.. ..- attr(*, "data_types")= chr [1:2] "N" "F"
..@ lines :List of 1
.. ..$ :Formal class 'Lines' [package "sp"] with 2 slots
.. .. .. ..@ Lines:List of 1
.. .. .. .. ..$ :Formal class 'Line' [package "sp"] with 1 slot
.. .. .. .. .. .. ..@ coords: num [1:2, 1:2] -61.5 -61.4 -34.4 -34.3
.. .. .. ..@ ID : chr "0"
..@ bbox : num [1:2, 1:2] -61.5 -34.4 -61.4 -34.3
.. ..- attr(*, "dimnames")=List of 2
.. .. ..$ : chr [1:2] "x" "y"
.. .. ..$ : chr [1:2] "min" "max"
..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
.. .. ..@ projargs: chr "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84
+towgs84=0,0,0"
Run of maptools::gzAzimuth:
gzAzimuth(from = l1.geo at lines[[1]]@Lines[[1]]@coords[1,],
+ to = l1.geo at lines[[1]]@Lines[[1]]@coords[2,])
[1] 46.52677
Run of geosphere::bearing:
bearing(p1 = l1.geo at lines[[1]]@Lines[[1]]@coords[1,],
+ p2 = l1.geo at lines[[1]]@Lines[[1]]@coords[2,])
[1] 46.65786
Angle given by Linear Directional Mean for l1:
47.53076
For l1.geo (in Lat/Lon) the ArcGIS' tool gives a warning that it it should
be projected but calculates this angle:
WARNING 000916: The input feature class does not appear to contain
projected data.
51.94714
I can provide further information or files but i think this can be enough
for now.
Thanks,
Regards
2016-05-25 4:32 GMT-03:00 Roger Bivand <Roger.Bivand at nhh.no>:
On Wed, 25 May 2016, Eduardo Diez wrote:
Dear everyone,
from ArcGIS.
I could find functions for performing a similar task but taking Origin ->
Destination points and only in Lat/Lon (geographic coordinates) namely:
geosphere::bearing and maptools::gzAzimuth. (Somehow they give different
results for the same set of points: around 0.05 degrees, maybe because of
the trigonometry)
Because of the nature of my work it has to be in UTM.
Interesting question. Could you please provide a test set of pairs of
coordinates with the fully qualified UTM projection (which datum in
particular, and if WGS84, which version - there are several), the output
from ArcGIS (which version), and the geographical coordinates you see in
ArcGIS (to remove the possibility that it is a projection issue on the R
side). They don't need to mean anything, only that they can show us what
the problem is in code and numeric terms. A link to an R file and an RDS
file would be helpful.
Roger
Perhaps someone skillfull in trigonometry can figure it out but i'm not the
case.
Is there a way to do this in R?
Thanks,
Eduardo
[[alternative HTML version deleted]]