Guest User

Untitled

a guest
Jul 17th, 2018
71
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 1.27 KB | None | 0 0
  1. library(ggmap)
  2. library(rgdal)
  3. library(ggplot2)
  4.  
  5. #####################################
  6.  
  7.  
  8. ## build data frame with triangulation data
  9.  
  10. t1 <- c(534325, 4858925, 338, 0955)
  11. t2 <- c(534383, 4859261, 290, 1010)
  12. t3 <- c(534386, 4859011, 301, 1015)
  13.  
  14. df <- as.data.frame(rbind(t1, t2, t3))
  15. colnames(df) <- c("Easting", "Northing", "Azimuth", "Time")
  16. df$Time <- as.character(df$Time)
  17.  
  18. ## plot coordinates with triangulation
  19.  
  20. ggplot(df, aes(x=Easting, y=Northing, label=Time)) +
  21. geom_point() +
  22. geom_text(nudge_x = 50) +
  23. geom_spoke(aes(angle = Azimuth, radius = -500)) +
  24. theme_bw()
  25.  
  26. ## convert utms to lat long
  27.  
  28. utmcoor <- SpatialPoints(cbind(df$Easting,df$Northing),
  29. proj4string=CRS("+proj=utm +zone=12"))
  30. lonlatcoor <- as.data.frame(spTransform(utmcoor,CRS("+proj=longlat")))
  31. colnames(lonlatcoor) <- c("lon", "lat")
  32.  
  33. df$lon <- lonlatcoor$lon
  34. df$lat <- lonlatcoor$lat
  35.  
  36. ## plot on base map
  37.  
  38. zoom <- 15
  39.  
  40. meanlon <- mean(df$lon)
  41. meanlat <- mean(df$lat)
  42.  
  43. basemap <- get_googlemap(center=c(lon=meanlon, lat=meanlat), zoom=zoom,
  44. maptype="hybrid")
  45.  
  46. ggmap(basemap) +
  47. geom_point(aes(x=lon, y=lat), colour="red", size=2, data=df) +
  48. geom_spoke(aes(x=lon, y=lat, angle = Azimuth), data=df, radius=-200) +
  49. theme_void()
  50.  
  51. Warning message:
  52. Removed 3 rows containing missing values (geom_segment).
Add Comment
Please, Sign In to add comment